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SUMMARY 


Motion planning and control for the joints of flexible manipulators has 
produced delightful results in the planning of joint motions to achieve specified tip 
motions. This work is performed principally by Dr. Dong-Soo Kwon, who has just 
successfully defended his Ph.D. Thesis. Dr. Kwon's work took the specified tip motion 
and used it as the input to the inverse dynamic model which calculates back to the 
necessary joint motions. Inverse dynamic models for flexible arms have been 
developed by other researchers, but our work has overcome significant handicaps of 
their research and revealed some basic properties of the inverse models of these 
systems. 

The flexible manipulator model of tip motion arising from joint actuation is 
nonminimum phase system. In linearized form this appears as a system with a zero of 
the transfer function in the right half plane. The inverse of this system thus has poles in 
the right half plane, a condition usually associated with instability for physical systems 
which must be causal. The alternative interpretation of right half plane poles is as an 
acausal system. An acausal system responds to an input before the input is experienced 
by the system. The difference in the unstable causal and the stable acausal is the region 
of convergence of the Laplace transform. While all physical systems, including the 
flexible arm, must be causal, the inverse dynamic system is not a physical system. It 
represents a calculation of what joint torques would need to be to achieve tip motion. 
An acausal inverse system tells us that we must begin actuation of the joint before the 
tip of the robot begins to move. In fact this motion must begin at a time of negative 
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infinity for the motion to be exactly as specified. An extremely good approximation for 
a flexible arm requires only a fraction of a second of hub motion before the tip of the 
arm begins to move. This advance action bends the link into a shape that at t=0 retains 
the tip at the original location. After t=0 the tip begins to move, and can in fact move 
with an almost arbitrary acceleration profile. To make the motion practical a smooth 
profile is chosen, which none the less effectively uses the maximum torque of the motor. 
The smoothness of the tip acceleration profile reduces the peak torques required and 
the order of the inverse dynamic model necessary to accurately command the physical 
system. Since Dr. Kwon has produced an extremely efficient linearized form of the 
inverse dynamics algorithm, it can be applied in real time to the tele-operated case, 
rather than requiring the complete rest to rest trajectory be precalculated with great 
computing expense. When an obstacle or the object we desire to contact is sensed, the 
inverse dynamic model plans the remaining trajectory and can complete the motion for 
the human operator. After contact the flexible states of the arm can also be prescribed 
to assume some static deflection. This creates a very natural way to apply force with a 
flexible arm which is superior to commanding joint motion alone, since joint motion is 
very small for a large tip force, even with a flexible manipulator. Dr. Kwon has 
experimentally demonstrated his approach for a single link device, and we plan in the 
coming grant period to apply it to our large two link arm. Mr. J.D. Huggins is 
developing a Ph.D. proposal to extend Dr. Kwon's work in that direction. Dr. Kwon 
will be working for Oak Ridge National Laboratories in the coming year on problems of 
long arms. 
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Coordinated motion of a tip manipulator and the large flexible manipulator 
carrying it enables the inertial forces of the tip manipulator to damp the vibrations of 
the flexible manipulator. This concept has been demonstrated in simulation and a 
major thrust of the current grant period has been to experimentally demonstrate it. The 
small arm we use, SAM (Small Articulated Manipulator), is electrically actuated in three 
degrees of freedom. We are only using two degrees of freedom in the initial studies of 
motion in a plane. Initial plans for our research used three Motorola 68000 
microprocessors operating in parallel to control SAM. The algorithms were 
programmed and initial motion of SAM was achieved, but the internally developed 
control hardware proved to be unreliable. Consequently the controller was moved to a 
commercially available AT bus, 80386 processor. This unfortunate delay has been 
compounded by mechanical problems as well, which are presently overcome. SAM is 
again under control under the new processor and the controller algorithms are being 
studied experimentally. Simulation results show that the inertial forces are effective in 
damping the flexible vibrations but cannot act with the control authority that joint 
actuation can. This approach has the advantage that higher bandwidth may be possible 
with the small arm than with the large arm actuators. Consequently higher modes can 
be effectively damped. 

Disparate cooperating manipulators differ in size and other characteristics. Mr. 
Jae Lew is studying new ways to take advantage of these potentially complimentary 
features. Under his scenario the large arm carries the payload and the small arm. The 
small arm, capable of precise position control, pushes on the environment and deforms 
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the large arm slightly to achieve the final payload placement. The large arm tries to 
reduce the force that the small arm must apply to it with a force control. In addition to 
the greater position resolution naturally available from a small arm, its location near the 
payload avoids the non-minimum phase problems of the noncolocated large arm joints. 
Non-colocation leads to stability problems when using high gains to achieve a 
responsive control. A master/ slave control law has been shown in simulation to 
successfully accomplish the task. The control used in the past year's work acted to 
decouple the large and small arm behavior. The initial simulations assumed that the 
force between the two arms could be measured, which may not be the most practical 
mechanical design. Mechanical design similar to our laboratory large/ small arm 
combination makes measurements at the tip of the small arm, where it contacts the 
environment, more practical. A new round of studies is underway with a control that 
does not rely on the decoupling control. We expect the simpler mechanical design from 
relocating the force transducers to justify the control modifications. 

Also explored in the past months is the effect of disturbances generated on the 
large flexible arm (RALF) as a subsequence of completing its task. An abrasive cutoff 
saw was mounted on RALF and the spectrum of excitation and the resulting motion 
was studied experimentally. 

Motion planning for robots in free fall has obvious implications for all articulated 
objects in space. The non-holonomic characteristics of these systems results in 
reorientation of the vehicle (robot) when the joints move through a series of positions 
even though the joint angles return to their original values. The bad news is that the 
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normal operation of the robot in free fall will cause the orientation of the robot to shift 
with respect to the object it works on. The good news is that desired adjustments in 
orientation in the robot can be achieved using such cyclical motions without the use of 
thrusters, based on the same principles that a cat uses to reorient itself when dropped 
upside down. This project thrust seeks to avoid the problems and capitalize on the 
advantages with proper planning of robot motion. Mr. Jonathan Cameron, on leave 
from the Jet Propulsion Laboratory, has completed his Ph.D. thesis proposal and is 
setting up his planning system. His initial approach is to use symbolic processing of the 
dynamics equations to establish suitable forms. For a given system a catalogue of 
motion types will be established. The index to this catalogue is crucial to using the 
results in a planning system. Motion of a given type will allow reorientation from one 
class of poses to another. The detailed path of motion to bring one from an arbitrary 
initial orientation to the desired final orientation will depend on refining the parameters 
of the motion and perhaps iterative simulation. Consequently, very efficient simulation 
equations are desirable. Robots on orbit are the most obvious application of the 
research results, but robots on earth that fall are also candidates. A robot that could use 
dynamic planning to change its landing orientation would be much less susceptible to 
damage on landing. 

In addition to the four research lines above supported by the NASA grant, 
additional work under way is of interest because of its potential future impact on 
NASA sponsored studies and because these students are exposed to and benefit from 
the NASA sponsored research in their own graduate education. 
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A camera system for end point measurement of the position of the robot relative to 
retroreflective landmarks is being integrated into our laboratory. This has the potential 
for providing information to be used directly in feedback control. The retroreflective 
targets selectively placed in the work space are easily and rapidly distinguished from 
the background even by our inexpensive camera system. We are currently limited to 
about 10 position measurements per second, fast enough to refine our commanded 
position values, but not to perform true end point position feedback. Klaus Obergfell is 
now writing up his Master's thesis on this topic. He is planning to continue at Georgia 
Tech on a Ph.D. program, in part because of the experience we were able to provide him 
in the environment of this NASA grant. The integrated camera and processing system 
is provided on a beta test basis by Dickerson Vision Technologies, Inc., a Georgia Tech 
spinoff company. 

In addition to the use of RALF as an autonomous arm, it is also commanded as a 
tele-operator. A small kinematic replica master arm has been created and the 
corresponding control algorithms created. Application of the improved flexible arm 
servo controls to the tele-operator case are under study. Dave Magee is midway 
through a master's thesis on this topic. He has just passed his Ph.D. Qualifying 
Examination and wishes to continue his research in similar areas. 

Research using high speed digital signal processing chips for control is also 
underway. Speeds of 20 MIPS are available, but to take full advantage of this speed 
some tailoring of the control algorithms is desirable. Lonnie Love is applying a 
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TMS320C25 processor to robotic exercise devices. Lonnie has also just passed his Ph.D. 
Qualifying Examination and is beginning his Ph.D. proposal on related areas. 

Finally, the effect of tapering the links of the arm are being considered by Doug 
Girvin. Doug is looking particularly at the nonminimum phase characteristics of the 
arm, i.e. the zeros of the transfer function from the joint torque to the tip translation. 
Since the traditional finite element and even assumed modes technique are 
approximations that are not verified in terms of their accuracy of zero location, this 
research is using an exact frequency domain technique. The exact solution is perfectly 
valid for the linear case only. Numerical evaluation of the exact analytical solution still 
introduces approximations, of course. By looking at various tapers we expect to be able 
to modify the minimum phase characteristic of the arm. 
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RESEARCH TOPIC: Control of a Flexible Bracing Manipulator 

RESEARCH ASSISTANT: Dong-Soo Kwon 

SHORTTERM OBJECTIVE: Integration of current research work to realize the 

bracing manipulator 

The tip position control of a flexible manipulator can be considered as the most 
general and complicated problem in the control of robotic manipulators, because it is 
the tracking control of the nonlinear, noncolocated, nonminimum-phase system with 
uncertainty. All research results about the flexible manipulator control were integrated 
to show a control scenario of a bracing manipulator. First, dynamic analysis of a 
flexible manipulator has been done for modeling. Especially, the noniminum-phase 
system characteristics have been studied in detail, and the relation between the 
nonminimum-phase system, and causal and anticausal systems has been analyzed. 
Second, from the dynamic model, the inverse dynamic equation is derived, and the 
time-domain inverse dynamic method were proposed for the calculation of the 
feedforward torque and the desired flexible coordinate trajectories. Third, a tracking 
controller is designed by combining the inverse dynamic feedforward control with the 
joint feedback control. The control scheme has been applied to the tip position control 
of a single-link flexible manipulator for zero and non-zero initial condition cases. 
Finally, the contact control scheme was added to the position tracking control. 
Therefore, a control scenario of a bracing manipulator is provided, and evaluated 
through simulation and experiment on a single-link flexible manipulator. 
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RESEARCH TOPIC: Control of a Small Working Robot on a Large Flexible 

Manipulator for Suppressing Vibrations 
RESEARCH ASSISTANT: Soo Han Lee 

SHORT TERM OBJECTIVE: Development of a robust control law for flexible robot 

and it's stability analysis 

1) Development of a robust control law for flexible robot 

In order to achieve the accurate trajectory control of a flexible robot, we need to 
adopt a nonlinear control law. We develop a robust control law for slow motions. The 
control law does not need larger velocity gains than position gains, which some 
researchers need to prove the stability of a rigid robot. Initial experimentation for SAM 
(Small Articulated Manipulator) show that the control law which use smaller velocity 
gains is more robust to signal noise than the control law which use larger velocity gains. 
The control law for the fast motion is strain rate feedback. 

2) Stability analysis 

When we adopt a composite control law for a two time scaled model, the rule of 
thumb is to separate the controller bandwidth of the two control laws at least factor of 
3. However, we may need a systematic approach for determining the controller gains 
which are related the bandwidth of the controller. Hence, we analyze the stability of 
the composite control law, the robust control for the slow motion and the strain rate 
feedback for the fast control. The stability analysis is done by using the quadratic 
Liapunov function. 
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3) Phase lag of input force 

The force generated by SAM effect the lower link flexible motion of RALF (Robotic 
Arm Large and Flexible) through upper link when the SAM is located at the tip of the 
upper link. One relationship between input force and the flexible motion of the root of 
upper link is nonminimum phase which some researchers have characterized the time 
delay. The other relationship is minimum phase because the force is parallel to upper 
link (no phase lag). These relationships are analyzed by assuming a time delayed 
model. 

We can control the flexible motion of links by relating the input force to the flexible 
signals which are sensed at the near tip of each link. The signals are contaminated by 
the time delayed input force. However, we can reduce the effect of the time delayed 
input force by giving a certain configuration to SAM. 

4) Present research 

The hardware and software are under development for the real time control of 
RALF and SAM. 
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RESEARCH TOPIC: 
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Control Strategy for Cooperating Disparate Manipulators 
RESEARCH ASSISTANT: Jae Y. Lew 

SHORT TERM OBJECTIVE: Non-Colocated Control of Disturbances of a Flexible Arm 
The long term objectives for this research use a small arm (SAM) to compliment a 
large arm (RALF). As a step on this path the disturbances generated by the tasks of the 
robot were experimentally studied during this period. In general, disturbances could 
arise from the motion of the payload (or another arm it carries), from planned or 
unplanned contact with the environment, or from processes at the tip of the arm 
(grinding, cutting, spraying, etc). The available equipment made it convenient to study 
disturbances generated by an abrasive cut off saw mounted on the tip of RALF and 
used to cut through rods and pipes in the range of 0.5 in. to 1.5 in. The abrasive process 
was modeled as a friction operation. 

The operation of the cut off saw was performed autonomously and under tele- 
operated control. The behavior was stable in both cases for a P-D joint control 
algorithm. The critical aspect of operation was low speed on initial contact. The angle 
of approach relative to the stiffness ellipsoid showed some significance, but since the 
ellipsoid was not greatly elongated leads to only minor differences for different 
approach directions. The broad band excitation during the cutting operation is capable 
of exciting arm natural frequencies to above the third natural frequency of the arm. The 
direction of cut changed the modes excited due to the change in coupling coefficients. 
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PUBLICATION 


1. Lew, J., Huggins, J. D., Magee, D. and Book, W. , "Experimental Investigations of 
the Effects of Cutting and on Chattering of a Flexible Manipulator." 1991 ASME 
Winter Annual Meeting, Atlanta, GA 30332. 



RESEARCH TOPIC: 


Dynamics of Multibody Systems in Free-fall 


RESEARCH ASSISTANT: Jonathan M. Cameron 

SHORT TERM OBJECTIVE: Complete and Defend Ph.D. Proposal 

PERIOD OF COVERAGE: January-March, 1990 (Winter Quarter) 

During this period, the researcher completed and presented his Ph.D. proposal. 
The researcher was only supported during winter quarter. After winter quarter, he was 
transferred to another project and was not supported by this contract. 

The Ph.D. dissertation proposal was titled, "Planning and Executing Motions for 
Multibody Systems in Free-Fall." It was presented and successfully defended on March 
22, 1991. The proposal is attached. Although an extensive literature search was 
reported in earlier reporting periods, more was done in preparation of the thesis 
proposal during this quarter. 

The researcher also supported activities in the Mechanical Research Building by 
assisting with computer system management required for several DEC VAX computers. 
PUBLICATION 

1. Cameron, Jonathan, "Planning and Executing Motions for Multibody Systems in 
Free-fall", Ph.D. Dissertation Proposal, George W. Woodruff School of Mechanical 
Engineering, Georgia Institute of Technology, Atlanta, GA 30332. 
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ABSTRACT 

Many applications of robotic and teleoperated manipulator 
arms require operation in contact and non-contact regimes. This 
paper deals with both regimes and the transition between them 
with special attention given to problems of flexibility in the links 
and drives. This is referred to as contact control. Inverse 
dynamics is used to plan the tip motion of the flexible link so that 
the free motion can stop very near the contact surface without 
collision due to overshoot. Contact must occur at a very low 
speed since the high frequency impact forces are too sudden to be 
affected by any feedback generated torques applied to a joint at 
the other end of the link. The effect of approach velocity and 
surface properties are discussed. Force tracking is implemented 
by commands to the deflection states of the link and the contact 
force. This enables a natural transition between tip position and 
tip force control that is not possible when the arm is treated as 
rigid. The effect of feedback gain, force trajectory, and desired 
final force level are of particular interest and arc studied. 
Experimental results are presented on a one link arm and the 
system performance in the overall contact task is analyzed. 
Extension to multi-link cases with potential applications are 
discussed. 

1. INTRODUCTION 

Control of the impact between the arnTi tip and its 
environment has been largely ignored in prior control research . 
Since contact forces may be high, it is important to understand the 
behavior well. This is especially true for operation with repeated 
contact such as when employing a bracing strategy or for material 
handling applications. Compared to the behavior before and after 
the impact, behavior during impact is more complex in several 
ways: Impacts give high frequency excitation to the arm and to 
the environment. Modeling these behaviors requires models for 
the distributed mass and compliance of the system. The contact 
force exhibits a highly nonlinear nature because it is unidirectional 
and dependent on relative dp and environment position. 
However, other complex phenomena important in free motion are 
not so important in contact control, for example the nonlinearities 
of the dynamic equations. Since impact intentionally involves low 
speeds and small position changes the dynamic nonline&rities are 
often insignificant. 

First, this paper presents a practical inverse dynamics 
approach for planning trajectories that will bring the tip of a 
flexible manipulator to a precise stop with no overshoot from a 


still or moving initial condition. Tip position as a function of time 
U uaed aa an input to the inverse dynamic equations, which 
generate the remaining states and joint torque input needed to 
achieve the desired tip position. 

After the planned trajectory is executed, contact with the 
surface is initiated. Planning of the tip position will allow the 
contact to be initiated from a point very near the surface. Since 
the deflection variables of the arm are modeled and controlled 
along with joint variables, force and tip position control are both 
naturally accommodated by the same model. A smooth time 
history of tip force is readily converted to a static displacement of 
the flexible modes of the links. 

The overall system performance involves time penalties, 
actuator constraints, and force constraints on the tip and the 
environment. Several characteristics of the application will favor 
the present more advanced approach. Applications with long arms 
and penalties on the structural size and mass clearly favor this 
advanced approach, especially when the surface contacted is prone 
to damage. Airplane washing, painting and depainting robots have 
been constructed and have these characteristics. Nuclear waste 
remediation inside large single shell storage tanks are another 
clear example, a a are required decontamination and 
decommissioning of nuclear facilities. Manufacturing and 
construction of large items such as aircraft, power generation 
equipment, bridges, and buildings also lend themselves to this 
technology. 

The organization of the paper will now be described. The 
inverse dynamics approach will be presented in summary form 
with reference to details in section 2. It forms the core of the first 
phase of contact control, motion control. The second phase is the 
impact phase, covered in section 3. Hie final phase of force 
tracking is then described, followed by a section presenting the 
experimental results on all phases. Section 5 collects the results 
from the various phases to show the effect on the overall task 
performance time. Section 6 presents the issues involved in 
extending the approach to a multi-link case. The final section on 
conclusions summaries the results. 

2. AN INVERSE DYNAMICS APPROACH WITH 
PRACTICAL COMPUTATIONAL EFFICIENCY 

To minimize the impact force and the task completion time, 
the position control should bring the tip near the contact wall 
precisely without overshoot. A practical and simple inverse 
dynamic method, which can make that motion possible, was 
proposed by Kwon and Book [4]. The inverse dynamic method 
generates the torque profile and all state trajectories that can make 



the tip of a flexible manipulator follow the desired trajectory. The 
trajectories are implemented with a tracking controller. The 
relative position of the surface is assumed to be poorly known 
when the motion begins. The initial condition depends on the 
state of the arm when a position sensor on the tip of the aim 
determines the distance to a surface in the environment. Upon 
determining the position of the environment the fast inverse 
dynamic algorithm allows the remainder of the arm’s trajectory to 
be planned on-the-fly , in time to complete the motion of the tip 
with precise position and no overshoot. This section describes the 
modeling of a flexible manipulator and the inverse dynamic 
method briefly. 
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2.1. Modeling 


A single link flexible manipulator having planar morions is 
described as shown in Figure 1. The rotating inertia of the servo 
motor, the tachometer, and the clamping hub are modeled as the 
hub inertia lb. The payload is modeled as the end mass Me and 
the rotational inertia Je. Though structural damping exists in the 
flexible link, it is ignored in modeling. 

To derive equations of motion of the manipulator, we 
describe the position of a point on the beam with virtual rigid 
body morion and flexible deflection using a Bemoulli-Euler beam 
model. The virtual rigid body motion is represented by the 
motion of a moving coordinate, which is attached to the beam. 
The flexible deflection is described by a finite series of assumed 
modes with respect to that moving reference frame. The rigid 
body morion can be defined differently according to the choice of 
the rigid body coordinate. Several authors [2,3] used the rigid 
body coordinate that is attached at the base hub or that passes 
through the center of mass of the beam. In this paper, the rigid 
body mode coordinate that passes through the end point of the 
beam is selected [1,4], and the mode functions of pin hub-inertia, 
pin end-mass boundary conditions are used to describe the 
displacement of the beam. Because of the rigid body coordinate 
selected, the end point position of the beam can be expressed by 
the rigid body mode variable alone. With this simple 
representation of the end point position, we can derive the inverse 
dynamics equation easily. 

By using Lagrange’s equations of motion, the dynamic 
equation of a flexible manipulator is obtained with generalized 
coordinates. 
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The mass, stiffness, damping, and input matrices are derived in 
these general forms. They are valid for different definitions of 
rigid body modes. The damping matrix represents the joint 
friction modeled as the viscous friction with the coefficient c 0 . 

The dynamic equation can be divided into a rigid body 
motion part and a flexible motion part as follows: 
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For a state space form, we obtain the following dynamic equation. 
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2.2. Inverse Dynamic Equations 

From the above the direct dynamic equations, we will derive 
the inverse dynamic equation. Eqn. (1) can be written in two 
parts. 
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From Eqn. (3), torque is expressed as 
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Substitution of above Eqn. (3a) into Eqn. (4) gives the following 
relations between the flexible coordinate and the rigid body 
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From Eqn. (4), the acceleration of flexible coordinate is expressed 
as follows: 


The inverse dynamic equation is obtained in this simple form. 
By the way, it seems at first to be impossible to integrate that 
equation because the inverse system matrix has positive real 
valued poles (which came from the positive zeros of the transfer 
function of the direct dynamic system) as weD as negative poles. 
However, if we relax the solution range of that equation to include 
noncausal solutions, we can obtain a unique stable solution by 
integrating this 1st order differential equation. 

The inverse dynamic system equation can be divided into the 
causal part and the anticausal part by using a similarity 
transformation as follows. 

[7]: Orthogonal transformation matrix 
- [T C .TJ{P C ,PJ T 


4, - 

. < 6 ) 

Substitute this Eqn. (6) into Eqn. (3), then we will get the 
following Eqn. 
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If Eqn. (5) and Eqn. (7) arc represented in a state space form, the 
inverse dynamic equations will be written in this form. 


where 4y} r , T * s basis includes the eigenvectors which 

have negative eigenvalues and T m is made of the eigenvectors of 
positive eigenvalues. 



Such a coordinate change decouples the inverse system into two 
subsystems of Eqn. (9). The new variable P c represents the 
coordinates of the causal system, and the P m represents the 
anticausal system. 

For a given end point trajectory of Figure 2, the causal part 
torque is obtained by integrating the causal part inverse dynamic 
equations forward from the initial time of the trajectory. Then, 
the anticausal system equations must be integrated backward from 
the final time of the trajectory. The total torque, which is the 
output of these equations, is obtained by adding the both part 
torques as shown in Figure 3. As additional outputs, the reference 
trajectories of the output of Figure 4 were obtained from the 



integration of the invene dynamic equation. 

3- CONTACT CONTROL 

Impact Force Impact is a high frequency phenomena dependent 
mainly on the local parameters at the point of impact, especially 
the tip velocity and the elastic, mass and damping properties of 
the tip and the environment. The delay inherent in the tip 
response any elastic link will defeat a strategy to reduce impact 
force by tip force feedback control of joint actuators located at the 
other end of the link. Only for an extremely soft environment 
will the response be fast enough. Low impact force and quick 
task completion requires accurate tip placement before impact and 
low velocity upon impact Another approach that deserves 
attention is to cushion the impact between aim and environment 
with a "shock absorber” system. This might be designed in 
conjunction with the improved control approach. 

To study impact, the flexible model of the link must be 
supplemented with a model of the surface compliance and 
damping. The dynamic model of the arm remains essentially the 
same as without contact, except the contact surface and the force 
sensor must be incorporated as additional stiffness and damping 
matrices. The tip position is readily available with the choice of 
coordinates used in the model. The tip is allowed to contact a 
surface modeled by a spring and damper. The tip is free to 
rebound from the surface as dictated by the combined dynamics, 
where upon the contact force drops to zero. 

Figure 6 is an overlay of the dp forces for various dp 
approach velocides. The initial impact results in rebound, loss of 
contact, and further impacts. Higher dp approach velocities mean 
higher impact forces but the change in velocity has minimal effect 
on the time between initial contact and the occurrence of the peak 
force which typically occurs on the first or second impact. Figure 
5 shows an expanded view of a single impact for various surface 
stiffnesses and displays the time scale of these impacts, which is 
of the order of .001 seconds. To reduce these impact forces with 
feedback force control, the control action would have to be taken 
and propagate to the tip in less than 1/2 the duration of the 
impact. This is clearly infeasible for the light weight arm, whose 
natural frequency is only 2.1 Hz. Even supposedly rigid arms 
will not have the speed of response needed to react in the required 
time because of the actuator limit. This leaves the designer with 
two alternatives: 

1) reduce the speed on impact, or 

2) modify the compliance at the point of contact. 

Either of these alternatives might be useful in particular situations. 

A "shock absorber” could be mounted on the tip to reduce 
the peak force. While the peak force has been reduced, the total 
momentum transferred to the environment would remain about the 
same, which may indicate problems for other parts of the system. 
If enough travel of the shock absorber is permitted, the arm might 
now be able to respond to the initial impact. The soft suspension 
allows the tip to be used as a proximity senior. Engineering this 
solution might limit the utility of the arm, however. 

The solution discussed in this paper is to reduce the speed on 
impact by carefully controlling the tip position. This requires that 
the planned motion trajectory bring the tip very near the surface, 
or else the low approach velocity will delay the progress of the 
task. The inverse dynamics approach allows this to be done for 
the flexible arm. The position of the surface can be updated as it 
comes within range of a proximity sensor and the trajectory can 


be modified en route. It is also possible to plan for a non-zero tip 
velocity at the end of the motion trajectory which becomes the 
approach velocity. 

force Tracking Since the deflection variables of the arm are 
modeled and controlled along with joint variables, force and tip 
position control are both naturally accommodated by the tame 
model. A smooth time history of tip force is readily converted to 
a static displacement of the flexible modes of the links. These 
displacements are commanded to achieve the desired tip force. 
The time allowed to reach the final desired force and the force 
tracking gains are important parameters in obtaining successful 
force control. 

The force trajectory used in this paper is s cubic polynomial 
fitted between zero and the desired final force. The force control 
mode is entered upon the initial contact. From the known desired 
force the necessary flexible state variables are computed, assuming 
in our case that the system is in static equilibrium. (Refined 
assumptions regarding the dynamics of the constrained system are, 
of course, possible. Referring to Eqn.2, one can see that the 
flexible variables ty are found to be 

», ■ jg 

qy-C'-'iStt 

where F* is the desired contact force, St, is the strain at X i , and 
<r is the submatrix of the output matrix C. 

As can be seen from Figure 6, the approach velocity has 
negligible effect on the force tracking and regulation. The desired 
force level also has little effect. The gains used in the force 
tracking algorithm do affect force tracking, however, as does the 
force trajectory as illustrated by shortening the time for achieving 
the desired final force level. 

The force control algorithm prescribes desired values for the 
contact force, the flexible states, qp and the joint variables. The 
tip contact force is fed back with a proportional control, and the 
flexible states and the joint variable is fed back with t PD control. 
Integral action could be added to eliminate the steady state error 
of imperfections such as friction. These imperfections are not 
included in the model, so the simulation shows no problems in this 
regard. This is not the case for the real system as observed later 
in the physical experiments. Increases in proportional gain are 
shown in Figure 7. Increased gains ultimately lead to instability. 

The force trajectory determines the abruptness of the changes 
that the flexible variables will be tracking. Figure 8 shows the tip 
forces resulting from shorter trajectories leading to the same final 
force value. The force trajectory is initiated on the first contact. 
By the time the arm has returned for the second contact the force 
trajectory may be half completed, as in the first case in Figure 8. 
This leads to very high force levels on the second contact. Even 
for the case of zero approach velocity, with the tip in zero force 
contact with the environment, the sudden force command can 
cause the tip to leave the surface due to a nonminimum phase 
reverve action. [2] 

4 . EXPERIMENTAL RESULTS 

The contact controller was applied to a single link 
experimental arm to verify the results. The parameters of the 
system are the same as for the simulated system. These 
parameters are shown in Table 1, which references the description 


of Figure 1. Figure 9 shows schematically the control system 
implemented. The measurements in the vector Y include joint 
position and velocity and strain at the base of the link and near its 
mid-point. The results of position and contact force control are 
represented in Figure 11 ,12,13- The desired, simulated, and 
experimentally measured contact force, tip position, and link base 
strain are displayed. Tip position, displayed in Figure 11, is 
calculated based on the joint position and the strain values, since 
tip position measurements are not directly available in this 
experiment. 


Table 1. Physical properties of the flexible manipulator 


Link 

El: stiffness ibfin* 

4120 

H: thickness in 


Ro A: unit length 
msss lbfcVin 3 

2.54E-4 

*0.14 

L: length in 

,7 

Tip miu 

Me: mass Ibm 

0.1 

Jerrot. inenis 

3.8E-5 H 

Hub 

Or. rot. inertia Ibjrtn 

0.14 || 


The experiments duplicate to acceptable accuracy the 
simulation results with the exception of the steady state error in 
the contact force of Figure 12. This steady state error is 
apparently due to the substantial static friction in the motor. The 
primary imperfection of the model used in simulation and 
controller design is the absence of friction. The friction is 
compensated for in the position control phase of the task. The 
friction model is not able to compensate for both the dynamic 
friction during motion and the static friction when the arm is 
almost stationary. Since there is a limit on the proportional gains 
for a stable response, an integral action is proposed to improve 
steady state error. 

Finally, the strain at the base, the desired, the simulated, and 
the experimental values, are shown in Figure 13 for all three 
phases,position control, approach and impact, and force trajectory 
control. The motion in this case is a 70 degree rotation of the 47 
inch long beam, that is an arc length of about 57 in. 

5. SYSTEM PERFORMANCE 

While a optimization of the overall performance has not yet 
been attempted, the trends permitted by the planning and control 
approach described here are clear. The complexity resulting from 
the additional sensors needed, the complexity of the control, and 
the need for reasonably accurate models of the arm are all reasons 
to use conventional approaches to contact control when they will 
suffice. 

The task time is broken into several parts: 

1) position tracking 

2) position regulating (tip position settles) 

3) velocity control (move toward impact) 

4) contact force tracking, and 

5) contact force regulating (force reading settles) 

The time for each part is shown in Figure 10. Position tracking 
is reduced by controlling the arm as a flexible system. Faster 
moves are then possible. Since the inverse dynamic approach 
results in zero overshoot, the position regulating pari is reduced 
to almost zero. Velocity control time is reduced since the inverse 
dynamic approach avoids overshoot. Contact force tracking time 


may be increased with the proposed method to bring about a 
reduction in the regulation or settling time. 

6. EXTENSION TO THE MULTI-LINK CASE 

The scale of the initial simulations and experiments justifiably 
lead one to question the applicability of the approach used here to 
multiple links and axes of motion. The most critical questions 
surrounding the technique are: 

1) Can the technique be extended to tip motions in three 
dimensions, and by what techniques? 

2) The single link case is well approximated by a linear 
model while an articulated multi-link arm has prominent 
dynamic nonlinearities. How limiting is the dependence on 
linear techniques? 

3) The number of states of an arm grows with the number of 
axes. A six degree of freedom robot would have at least 
three axes and links to be considered for positioning a contact 
point. How does the computational complexity scale with the 
number of degrees of freedom. 

The approach, even in a single degree of freedom bears 
promise for high performance servo applications. In addition, 
there it no conceptual difficulty in extending the techniques from 
one dimension to three dimensions. The mode shapes used to 
model the flexible beam should continue to be of the same 
boundary conditions to enable position to be found as a function 
of the rigid variable alone. Other issues, such as a flexible 
equivalent of hybrid control appear possible, but the detailed 
analysis has not been carried out. 

The dynamics of a multi-link arm involve nonlinear terms not 
found in the single link case. These nonlinearities are the second 
order velocity terms (Coriolis and centrifugal forces), the 
variations in the inertia matrix (joint angle dependent), and the 
variations in the gravity forces (joint angle dependent). The 
contact control problem arises under the relatively low velocity 
situation near contact. The square of the small velocities will 
therefore be even smaller. Since the velocities are small, the 
position changes will also be small. A long slender arm 
supplements these arguments by moving small amounts in terms 
of joint angles for a large tip motion. Our large two degree of 
freedom experimental arm, for example, has been verified to show 
little nonlinear behavior even at maximum joint velocity. [5] 

The final issue in the extension is the calculation time. As a 
last resort it is possible to calculate the motion trajectory off-line 
before hand. The current trajectory is calculated before the 
motion on a MicroVAX Q computer, roughly a 1 MIP machine, 
in 1.01 seconds for 57 inch long tip trajectory (1.066 sec travel 
time). RISC processors and DSP chips with 20 MIP performance 
are commonplace and inexpensive. When progressing from one 
link to three links, how does the computational complexity scale? 
For the one link case of this paper, six state* are found in the 
flexible model. The cost effective method of inverse dynamics 
reduces this to four because the nominal (rigid body) dynamics are 
specified by the desired tip position. The remaining fourth order 
system is separated into two second order system by the causal 
and anti-causal property of the solution. The doubling of states 
causes twice as many state equations, each twice as complex. 
Hence the computational complexity grows as the square of the 
degrees of freedom. The increase by a factor of 9 in 
computational complexity is more than offset by the possible 
increase by a factor of 20 to 40 in computational power readily 



available. It is also possible that the algorithms themselves could 
be further streamlined. 

7. CONCLUSIONS. 

An effective scheme for contact control of a flexible link arm 
has been derived, analyzed and experimentally verified. It 
provides fast tip motion with no overshoot, impact with modest 
impact forces, and force tracking up to a desired constant force 
level. There are several avenues to improve and extend the 
method. The primary one discussed here is the extension of the 
technique to multiple degrees of freedom. It appears that this 
extension is feasible with current technology. 

Applications of the approach to aircraft servicing, nuclear 
waste remediation, and manufacturing of large parts seems to be 
among the most suitable applications at this time, since long, light 
weight arms are involved and there is a need for repeated contact 
with the work piece. 
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Figure 2. Desired end point trajectory 
a) Acceleration, b)Vek>city, c)Position 



Figure 3. Calculation of torque with the 
inverse dynamic method 
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Figure 4. Trajectory generation by 
using the inverse dynamic method: a) 
Joint angle, End point position, b) 
Strain at base and mid-point 


Figure 1. A single link flexible 
manipulator 





Figure 5. Impact force* for different contact 
surfaces 


Figure 8 Contact force for desired force 
tnjecories with variou application time 



Figure 6 Contact force for various tip approach 
velocities 

Figure 9. The control scheme of the experiment 



Figure 7 Contact force for difFemet force 
feedback gains 


Motion Planning 






Figure 10. Comparison of task completion times 
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Figure 11 The experimental results of the tip 
position 



Figure 12 The experimental result of contact 
force 
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Figure 13 The experimental result of the strain 
at the base 



N91-30528 

An Inverse Dynamic Tracking Control for 
Bracing a Flexible Manipulator 


SUMMARY 

A manipulator can improve its positioning accuracy in a manner similar to the 
human who braces his wrist for a more accurate finger motion. Reduction of positioning 
uncertainty by bracing can increase practical applications of flexible large-workspace 
manipulators. The goal of this research is to realize the bracing manipulator, which is 
flexible. In order to brace the manipulator against a surface of the environment, a tip 
position tracking control becomes essential to avoid large impacts on contact, and to 
control the force after contact. If we consider the characteristics of the flexible 
manipulator, the control task may be characterized as a tracking control of a nonlinear, 
noncolocated, nonminimum-phase system with uncertainties. 

First, the author introduces the staged positioning concept and provides a 
conceptual background of the bracing strategy. Second, in order to understand the 
dynamics and determine the control strategy, a single-link flexible manipulator is 
modeled by the assumed-modes method, and the validity of the model is verified with 
experimental results. Several techniques to obtain more accurate mode functions are 
discussed and compared. Third, a time-domain inverse dynamic method is proposed to 
cancel out the nonminimum-phase positive zeros, which cannot be canceled by feedback 
controls. The trajectory is also planned to use the maximum capacity of an actuator and 
to avoid exciting flexible vibrations. Fourth, the tracking controller is designed to 
control the end-point of a flexible manipulator without any overshoot or residual 



xiii 

vibrations. Its perfect tracking or asymptotic tracking performance is analyzed and 
discussed for zero and non-zero initial condition cases. Finally, impact phenomena are 
investigated, and the contact force control of flexible bracing manipulator is presented 
to set up a control scenario of a bracing manipulator. The proposed control scheme is 
implemented on a single-link flexible manipulator, and analyzed and evaluated through 
simulations and experiments. 
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Tracking Control of a Nonminimum-Phase Flexible Manipulator 
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Abstract 

Perfect tracking or asymptotic tracking of the tip position of a flexible manipulator , which 
is a typical nonminimum-phase system, has been very difficult because of the positive zeros of 
the nonminimum-phase system and the lack of desired trajectories of flexible coordinates. This 
paper presents a tracking control scheme combined with the inverse dynamic feedforward 
control. The inverse dynamic method calculates the feedforward torque that cancels poles and 
zeros of the nonminimum-phase system. It also generates the desired flexible coordinate values, 
which match equivalent to the tip position trajectory dynamically. The feedback loop achieves 
tracking capability with the calculated desired flexible coordinate trajectories. The control 
scheme has been applied to the tip position control of a single-link flexible manipulator for zero 
and non-zero initial condition cases. The non-zero initial conditions of the system states are 
divided into three components of rigid body, causal and anticausal parts. The rigid body 
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component is used for the desired tip trajectory planning and the other components of these are 
used separately for the calculation of the feedforward torque of causal and anticausal parts. 
Through simulation and experiment, we explore the effectiveness and limitation of this method 
for moving non-zero initial condition cases. 

1. Introduction 

More and more today, robot manipulators are being used in the manufacturing processes, 
such as welding, material handling, and assembly. An increasing number of applications 
demand the capability of larger work-space, more accurate positioning, faster motion, less power 
consumption, and greater payload. The lightweight, large work-space manipulator has intrigued 
the engineers who have to deal with large objects, such as the space station assembly, large 
structure welding, airplane assembly, and nuclear waste handling. However, the large 
lightweight structure results in inherent flexibility, which induces complex dynamics in a 
manipulator. 

The dynamic model is infinite dimensional due to infinite flexible modes. It is nonlinear, 
and always has uncertainty due to the unmodeled dynamics. The tip position control with the 
joint torque is usually nonminimum-phase due to flexibility, and noncolocated because of 
different sensor/actuator location. In spite of the above characteristics of a flexible manipulator, 
the end-effector position should accurately track the desired path for some practical applications. 
Therefore, our control objective becomes the tracking control of the nonminimum-phase, 
noncolocated, and nonlinear system with uncertainties. This may be the most general and 
complicated problem in the control of robotic manipulators. 


t 
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To design a proper control scheme for a flexible manipulator, we have to decide what 
to overcome, what to avoid, or what to compromise. Most recently proposed ideas for the 
control of a flexible manipulator can be classified into following categories. Each method seems 
to have its advantages and penalties. 

a) Use many sensors and actuators 

If the number of actuators and sensors are equal to the number of controlled modes, the flexible 
structure can be controlled easily without a spillover problem [14]. Crawley [8] showed the 
possibility by using a distributed array of piezoelectric devices for actuation on a flexible 
structure. This may be the direct approach to the multiple mode system, but the penalty is the 
hard ware cost. 

b) Change the structure of the manipulator 

Alberts, et al [1] used a thin film of viscoelastic material between the structure surface and an 
elastic constraining layer to obtain passive damping. Park and Asada [15] tried to change the 
nonminimum-phase flexible manipulator system to a minimum-phase system by relocating the 
torque actuation point. This approach may be limited because the modification of the structure 
is not always possible. 

c) Use advanced feedback control schemes 

Hastings and Book [9] used optimal control with the strain feedback. Cannon and Schmitz [6] 
demonstrated the end-position control with the tip position and the joint rate feedback. Kotnik, 
et al [11] presented the results of the end-point acceleration feedback. Wang and Vidyasagar 
[19] applied the stable factorization approach to obtain the optimal step response. Singular 
perturbation methods were tried by Siciliano and Book [16], Khorrami and Ozguner [10], and 
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others. Various adaptive techniques have been proposed by Yuh [21], Yuan and Book [20], 
Cetinkunt and Wu [7], and other researchers. Bartolini and Ferrar [3] tried to extend the 
adaptive pole assignment control to nonminimum-phase system. Even though many control 
schemes have been proposed to control the flexible vibration actively, most control schemes will 
allow bounded tracking errors because the feedback control cannot cancel the nonminimum-phase 
zeros that have positive real values. 

d) Allow noncausality 

Bayo and Paden [3,4] developed the frequency-domain inverse dynamic method to calculate the 
feedforward torque profile for the end-point tracking. Kwon and Book [12] presented a time- 
domain inverse dynamic method, which generates desired flexible coordinates trajectories as well 
as the torque profile. These methods cancels the nonminimum-phase zeros, but the penalty is 
the loss of causality. 

e) Modify the input to not excite the vibration 

Seering and his co-workers [17,18] showed the residual vibration can be significantly reduced 
by shaping input signal with a series of impulses. The input shaping technique is effective for 
the simple mode system, but it is very sensitive to the damping and resonant frequency accuracy 
of the dynamic model. 

Among the above control methods, since the inverse dynamic method can truly 
compensate the nonminimum-phase positive zeros, and can achieve the tracking performance, 
the asymptotic tracking controller is designed using the time-domain inverse dynamic method 
in this paper, and implemented on a single-link flexible manipulator for zero and non-zero initial 
condition cases. This paper shows that the inverse system can be divided into causal and 
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anticausal parts, and each system can be integrated separately with the initial conditions of the 
system states that are divided into rigid body, causal and anticausal parts. The rigid body 
component is used for the desired tip trajectory planning and the other two components are used 
separately for the calculation of the feedforward torque of causal and anticausal parts. 

In the following sections, first, a single-link flexible manipulator is described and 
modeled using the assumed mode method. Second, the time-domain inverse dynamic equation 
is derived from the dynamic equation. Next, a tracking controller for the nonminimum-phase 
system is designed combining the inverse dynamic feedforward control and the joint feedback 
control. The tracking convergence characteristics are discussed in the Laplace domain, and the 
simulation and experimental results are presented for various initial condition (I.C.) cases 
including zero I.C. and non-zero I.C. 


2. Modeling 

A single link flexible manipulator having planar motions is described as shown in Figure 


1. The rotating inertia of the servo motor, 



the tachometer, and the clamping hub are modeled 
as a single hub inertia Ih. The payload is 
modeled as an end mass Me and a rotational 
inertia Je. Though structural damping exists in 
the flexible link, it is ignored in modeling. 


Figure 1. A single link flexible manipulator 



6 


Table 1. Physical properties of the flexible manipulator 


Link 

El: stiffness (lb f in 2 ) 

4120 

H: thickness (in) 

3/16 

Ro A: unit length mass 
(lb f s 2 /in 2 ) 

2.54E-4 

*0.14 

L: length (in) 

47 

Tip mass 

M e : mass (lbm) 

0.75 

J e :rot. inertia (lbfS^n) 

3.8E-5 

Hub 

I h : rot. inertia (lbfS^n) 

0.14 


To derive equations of motion of the manipulator, we describe the position of a point on 
the beam with virtual rigid body motion and elastic deflection using a Bemoulli-Euler beam 
model. The rigid body mode coordinate that passes through the end point of the beam is 
selected [2,12], and the mode functions of pined hub-inertia, pined end-mass boundary 
conditions are used to describe the flexible displacement of the beam. Because of the selection 
of rigid body coordinate, the end point position of the beam can be expressed by the rigid body 
mode variable alone. With this simple representation of the end point position, we can derive 
the inverse dynamics equation easily. 

By using Lagrange’s equations of motion, the dynamic equations of a flexible manipulator 
are obtained with generalized coordinates. The detailed equations can be found in the reference 
[ 12 ]. 


[M]q + [D]q * [K]q = [JB] t 


The dynamic equations can be divided into a rigid body motion part and a flexible motion part 


as follows: 
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where q r - q 0 ; rigid body coord.. 
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<tf s Y 


; flexible mode coord. 


( 1 ) 


For a state space form, we obtain the following dynamic equation. 


X - 

0 I 


0 ' 

M~ l K M' l D 

X + 

M'B 



Y » [C]x + [F]t (2) 

where X - [q r ,q f q r ,q f ) T 


3. Inverse Dynamic Equations 

From the above the direct dynamic equations, we will derive the inverse dynamic 
equation, which represents the relation between the desired acceleration of the rigid mode 
(equivalent to the tip acceleration) as input, and the torque as output. Eqn. (1) can be written 
in two parts. 


[M rr ]q r + [MJfy [D n ]q r + [DJq f -\Bfr 
[M/q r +[MJqy[DJ%+[DJq f +[KJq r [B^ 


(3) 

(4) 


From Eqn. (3), torque is expressed as 
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t = [5,]' 1 { [M rr ]q r * [MJq f + [DJq r * [DJq f ) (3a) 

Substitution of Eqn. (3a) above into Eqn. (4) gives the following relations between the flexible 
coordinate q f and the rigid body coordinate q r . 

[M^q r [D^q r [K^q f = [BJq+lBJ^ 

where [M,] . [IMJ-lBfiBX'WJ] 

[£>,] - ([D J ]-[B / ][S r l- , [D^| (5) 

TO - [*y 

[B u ] - I[S^[SJ' l [OJ-[D^ I '| 

[fij - 

From Eqn. (4), the acceleration of the flexible coordinates is expressed as follows: 

q f = -\M^[M/q r [M^[D/q r 

( 6 ) 

Substitute Eqn. (6) into Eqn. (3), then we will get the following Eqn. 


t - [CJ q f + [C a ] q f ♦ [F l7 ] + [FJ q r 


(7) 
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where 

If Eqn. (5) and Eqn. (7) are represented in a state space form, the inverse dynamic equations 
will be written in this form. 


Let X i = {q f ,q f } T , q ir = te r ,4 r } 

0 0 




o 

•^1 


m:'d, 


l 


m; 1 b u M; l B a 


Wir 


1 = 


*1 = W + iw» 

T = [C.]X ( ■+* [Ftfi, 


( 8 ) 


Even though the inverse dynamic equation is obtained in this simple form, it seems at 
first to be impossible to integrate that equation because the inverse system matrix has positive 
real valued poles (which came from the positive zeros of the nonminimum-phase direct dynamic 
system) as well as negative poles. However, if we relax the solution range of that equation to 
include noncausal solutions, we can obtain a unique stable solution by integrating this 1st order 
differential equation. 

The inverse dynamic system equation can be divided into the causal part and the 


10 


anticausal part by using a similarity transformation as follows. 

[7]: Orthogonal transformation matrix 

x, - m*, 

- [T e ,TJ{P c >PJ T 

where X.- {q f , q f ) T , T c 's basis includes the eigenvectors which have negative eigenvalues and 
T m is made of the eigenvectors of positive eigenvalues. 




Such a coordinate change decouples the inverse system into two subsystems of Eqn. (9). The 
new variable P c represents the coordinates of the causal system, and the P ac represents the 
anticausal system. 

For a given end-point trajectory of Figure 2, the causal part torque is obtained by 
integrating the causal part inverse dynamic equations forward from the initial time of the 
trajectory. Then, the anticausal system equations must be integrated backward from the final 
time of the trajectory. The total torque, which is the output of these equations, is obtained by 
adding the both torques as shown in Figure 3. As additional outputs, the flexible coordinate 
trajectories, which give the joint angle and strain trajectories, were obtained from the integration 
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of the inverse dynamic equation. 



Figure 2. Desired end point trajectory Figure 3. Calculation of torque with the 

a) Acceleration, b)Velocity, c)Position inverse dynamic method 


4. Analysis of the System State Space 

It is very useful to analyze the relations between the dynamic system states, and causal 
and anticausal subsystems’ states to construct additional desired trajectories of flexible 
coordinates. We can obtain the reference trajectory of the flexible mode coordinates of the 
direct dynamic system from a given rigid mode trajectory and the calculated state trajectories 
of the causal and anticausal systems. 

As we can expect in Eqn. (2), (8), and (9), the space of the full state vector X of the 
direct dynamic system can be divided into three subspaces: the rigid body coordinate subspace 
q^, the causal part flexible coordinate subspace P c , and the anticausal part flexible coordinate 
subspace P M . These subspaces are linearly independent and orthogonal to one another. The 
relations of these spaces are illustrated in Figure 4, and described by the following Eqn. (10) 
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when only two flexible modes are considered. 


where X={q r , q flt q fl , q r , q fl , q^} 7 
Qir~ 9/2» £/z) 


1 o' 
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* ir 

0 0 0 0 

0 0 


0 0 10 

p 0. 


pool 


X = H r q ir + HjiTP') 

X = + /// C P C + 


( 10 ) 



Figure 4. Dimensional analysis of state 
variables of flexible manipulator dynamic 
equations 


From the given end point trajectory, the rigid 
body coordinate trajectory q ir is obtained. The 
flexible coordinate trajectories of P c and P ac are 
then calculated from the integration of Eqn. (9). 
Thus, the trajectories of the complete state vector 
X can be obtained by using Eqn. (10). These 
trajectory values can be used as reference 
commands for feedback tracking control. 

The generation of complete state 
trajectories is an advantage of this inverse 


dynamic method. We get the strain trajectories as well as the joint trajectory. We no longer 
have to give the flexible manipulator reference commands to follow the trajectory like a rigid 
manipulator. 
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5. Tracking Controller Design 

Let’s consider a linear nonminimum-phase model (the output is the tip position), and a 
linear minimum-phase model (the output is the joint angle) of a flexible manipulator in the form 
of 


. nonminimum-phase 
Z t (s) 


X t (s) = 


Pis) 


t is) , 


. minimum -phase 

Z ,« 


0 is) - 


Pis) 


t(j) 


where X e : Tip position 

0 : Joint angle 

x : Joint torque input 


(ID 


Pis) = s 2 n is ±jp t ) 

i-1 

z,(s) = n (s± z<1 ) 

i* 1 

z B is) = n is±jz Bi ) 

i-l 

The nonminimum-phase system has imaginary poles p; and real positive and negative zeros z, 
compared to the minimum-phase system, which has the same imaginary poles and imaginary 
zeros z ei . 

The control objective is to make the output X e (t) (the tip position) follow the desired 
time-varying trajectory X ed (t) applying the joint torque r(t). In Figure 5, the input torque r will 
consist of the inverse dynamic feedforward torque r d , the feedback torque r b driven by tracking 
errors, and the disturbance or friction torque w. 

x = x d +x b + w 


( 12 ) 
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Figure 5. Tracking control of a nonminimum-phase flexible manipulator 

The feedforward torque r d is obtained using the inverse dynamic method that is equivalent to the 
two-sided Laplace transform with the region of convergence along the imaginary axis to get the 
stable solution in the time-domain. The feedback torque r b is generated by the error signals 
between the desired joint and strain trajectories, and the measured joint angle and strains. In 
experiment, the joint angle is measured directly without being estimated from the end-point X e 
measurement. Therefore, in Figure 5, the actual feedback loop is the dot line (...) rather than 
the solid line ( — ). Using the colocated feedback control with the joint trajectory, we could avoid 
the noncolocated control problem of the direct the tip position feedback. This advantage is 
obtained from generation of the desired flexible coordinate trajectories, which gives the joint 
trajectories and- the strain trajectories that match the tip trajectory dynamically. 
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T * = 


P'(s) 


XJs) 


Z' t (s) 

T b = m(Q/s)-6(s)) 
Ze( 5 ) 


where 0 /.s) = 


ZfjLs) 


A/x Z 9 (s) /s Z' 6 (s) 
Q(s) = ——x(s) or -j- 

p ( s ) Z (s) 


X t (s) 


The ’ notation of Eqn. (13) means the estimated model of the real system. 
From Eqn. (11), (12), (13), the output X^s) will be in the form of 


(13) 


z, W[ Z e (P f + KZ'q) 
P + KZ 6 z' t (P + KZJ 


v ed 


(14) 


If the feedback gain of K(s) becomes very large, the disturbance effect will be negligible, so that 
the tracking performance will be depend on the accuracy of the zeros of the model as follows. 


zX 

zX 


'■td 


(15) 


From Eqn. (14), the tip position tracking error dynamics will be as follows: 


KZ a z P’Z Y Zo 

(1 + -)E = — W+ {( i - 1) + - (Z t — -Z^)} x ed 

p p p Z > P Z' * 


(16) 


where E = X, - X. 


If the disturbance effect is negligible, and the dynamic model is good enough to cancel the 
system dynamics, the error dynamics will converge to zero as shown in Eqn. (17). 
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(1 + -) £ = 0 If W~ 0 

P 

P' = P (17) 

K - z t 

Then, the tracking dynamics are determined entirely by the joint feedback loop characteristics, 
which is a stable colocated control loop. Therefore, the control scheme will guarantee perfect 
tracking if the initial tracking errors are zero, and asymptotic tracking convergence for non-zero 
initial tracking errors, as long as the model is accurate. Even if the model is not exact, the 
tracking error will be bounded because the positive poles of the inverse system give a stable 
torque profile using the two-sided Laplace transform. 

6. Perfect Tracking with Zero Initial Conditions 

The above combined control scheme of the inverse dynamic feedforward control and the 
feedback control was implemented for a rest-to-rest motion of the experimental manipulator. 
For a real time control, a Micro Vax II was used with 12 bit A/D and D/A boards. The off line 
calculation of a trajectory and a torque profile was also performed in the Micro Vax. 

By applying the precalculated torque, compensating the joint friction, and using the 
feedback of the tracking error at the joint, the excellent result of Figure 5.5 was obtained. The 
flexible manipulator could stop without any overshoot or any residual vibration after it moved 
40 inches (48.76 degrees) within less than 0.8 second. In the strain signal, there exists a rough 
jerk that could be eliminated by using a smoother acceleration profile. Unfortunately, since an 
end-point position sensor was not available for the system, the end point position couldn’t be 
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measured directly. However, the end point tracking performance can be estimated from the joint 
tracking result and the strain tracking result. If the joint doesn’t have any overshoot or vibration 
and the strain doesn’t show any residual vibration, the end point can be presumed to stop without 
any overshoot or vibration. 

The experimental results show that a simple joint feedback PD controller provides 
excellent tracking if it is combined with the inverse dynamic feedforward control, and if the joint 
trajectory are generated from the desired end-point trajectory considered flexible dynamics. Full 
state feedback has been shown to converge even faster, but the increased difficulty in 
implementation seems unwarranted. 



Strain 



Figure 6. Experimental results of the combined control of the inverse dynamic feedforward 
control and the joint tracking feedforward control: a) Joint angle, b) Strain at the base 

7. Tracking Control for Non-zero Initial Condition Cases 

The previous section showed good tracking results for rest-to-rest motion, which gives 
zero initial and final flexible coordinate values. How, then, can this method can be applied to 
non-zero initial condition (moving) cases? In the control scenario of Book’s bracing manipulator 
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[13], the large flexible manipulator will be commanded to move in free space using a simple 
feedback controller with rough information of the target surface position until it senses the target 
bracing surface. If it measures more accurately the distance from the tip position to the bracing 
surface during motion, the controller will modify or create the desired tip trajectory 
accommodating the moving non-zero initial conditions, and calculate the required torque and the 
desired trajectories of all states with respect to the new tip trajectory. Then, it will follow the 
trajectory and will stop just in front of the surface or approach the surface with a suitable slow 
speed, and will contact to the surface without a large impact. Therefore, we need to be able to 
do trajectory planning and tracking control for general moving initial conditions. 

7.1 Analysis of Initial Conditions of the States and Trajectory Planning 

As explained in Figure 4 of section 4, the system states X of a flexible manipulator’s 
dynamics can be divided into a rigid body component q^, a causal component P c , and an 
anticausal component P, c . Naturally, the initial system states Xo can be divided into q^, P c0 , and 

P»cO- 

*0 = HflirO + ^f^c^cO + 

First of all, let’s assume the initial condition (I.C.) X 0 and the desired final state X f is 
known or measured exactly. The final condition X f also will be divided into three components 
q^, P cf , and P*. f . From the initial and final values of the rigid body coordinate, we are 
generating the desired end point trajectories including the desired acceleration, velocity, and 
position profiles considering the minimum traveling time and the high frequency content. The 
desired trajectories of Figure 7 are generated with the following example data: the end point 
initial acceleration X^q = -10 in/s 2 , the initial velocity X^o = 60 in/s, the initial position X^ 
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= 16.4 in, and the desired final position X ef = 57.4 in, X ct0 = 0, X^o = 0. The final condition 
doesn’t have to be stationary: it can be an arbitrary moving condition. Now we have planned 
the desired tip trajectories of a flexible manipulator. If the tracking controller is working 
properly, the end point will follow the trajectory and will stop without any overshoot or 
vibrations. 



Figure 7. The desired end-point trajectories for non-zero I.C. 


7.2. Inverse Dynamic Torque Calculation 

Since the initial and final conditions, and the rigid body coordinate trajectory (equivalent 
to the end-point trajectory) are given, the torque profile can be calculated using the same inverse 
dynamic method for the zero I.C. case. The causal torque is calculated by integrating the causal 
part equation with the initial condition P c0 , and the anticausal torque is calculated by integrating 
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Inverse Dynamic Torque 
(Non-zero I.C.) 



Figure 95.7 Inverse dynamic torque for non-zero I.C. 

the anticausal part equation backwards with the final condition P« f as shown in Figure 8. That 
means the anticausal part of the I.C. P, c0 cannot be considered in the inverse dynamic torque 
calculation. This is a limitation of the inverse dynamic method for non-zero initial condition 
case with the current fixed acceleration profile of the end-point. If the acceleration has more 
degree of freedom with extra parameters, it can satisfy the anticausal part initial condition by 
manipulating the auxiliary parameters. However, this approach of changing the acceleration 
profile structure is not so practical if implemented. In this research, we maintain the 
acceleration profile of Figure 2 adjusting the acceleration and deceleration time ratio without 
trying to match up to the anticausal part I.C. Even though the initial conditions of anticausal 
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part flexible coordinates cannot be considered in feedforward inverse dynamic control, the 
tracking error can be made to converge to zero by the asymptotic tracking capability of the 
feedback loop. The following sections examine several mismatched cases of initial tracking 
errors. 

7.3 Simulation Studies of Several Initial Tracking Errors 
Case 1 : Perfect knowledge of I.C. (no initial tracking errors) 

In order to get zero initial tracking errors, we have to know and use the exact initial 
values of the rigid coordinate and causal part flexible coordinate for trajectory planning and 
torque calculation, and the initial values of the anticausal part flexible coordinates naturally must 
be the same as the calculated anticausal flexible coordinates at t=0, which were obtained from 
the backward integration. Then, perfect tracking can be expected as shown in Figure 5.8. The 
simulation results perfectly match the desired trajectories of the end-point, the joint angle, and 
strains. 

Case 2 : Initial errors in flexible coordinates 

The assumption of the case 1 is that the initial values of the anticausal flexible 
coordinates naturally meet with the calculated values of the backward integration. If we relax 
this unrealistic assumption or allow the measurements of the strain to have errors, we will 
always have initial tracking errors of flexible coordinates. As the anticausal part torque near 
t=0 cannot cancel the dynamics of the system, but it is good enough to cancel the dynamics at 
t=t f , Figure 5.9 shows the asymptotically converging results. Though the strain signals have 
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oscillating tracking errors initially, they converge to the desired trajectories, and produce no 
overshoot or residual vibrations at the final position. 

Case 3 : Initial errors in rigid body coordinates 

In practical situations, it can be assumed generally that the measurement of the tip 
acceleration, velocity and position may have noise or errors. In this case, the desired end point 
trajectory is generated using the mismatched initial conditions. The rigid body coordinate 
trajectory (equivalent to the end-point trajectory) based on the mismatched I.C., will affect the 
causal part and the anticausal part torque calculation. Therefore, the feedforward torque cannot 
achieve the desired tracking characteristics, but the feedback control will compensate the 
tracking errors asymptotically. As shown in Figure 11, the tip trajectory converges slowly to 
the desired trajectory. This example was simulated with a wrong estimation of the initial tip 
velocity 60 in/s for the actual initial tip velocity 100 in/s. Intentionally, large initial errors were 
used for clear illustration of the tracking convergence. 

With this case, we can expect the estimation of the initial value of the rigid body 
coordinate to be very significant for the inverse dynamic method, because it is used for 
trajectory planning and it affects the torque calculation directly. 
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Figure 10. Case 2: Asymptotic tracking with initial errors of flexible coordinates 
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Figure 11. Case 3: Asymptotic tracking with initial errors of rigid body coordinates 
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7.4 Practical trajectory planning method 

The previous simulation analysis is based on the assumption that any calculation of the 
trajectory planning and the torque and desired trajectory calculation can be performed within 
next sampling time just after measuring the initial conditions. In order to use this control 
method for practical application, we had better consider the fact that the calculation time is much 
longer than one sampling period of the controller. Moreover, the sudden application of the 
torque is not desirable because it contains high frequency enough to excite the unmodeled high 
system natural frequency. 

As a practical method, when the actual moving initial conditions are measured at a certain 
instance, we estimate the states after noncausal time t=to* of Figure 12, and calculate the end 
point trajectory from the position at t=to’ to the final position at t=tf. With this end point 
trajectory, the inverse dynamic torque and the desired trajectory has been calculated. To avoid 
the sudden step torque, the noncausal part torque between t=t 0 and t=to* is applied as the 


Inverse Dynamic Torque 


(Non-zero I.C.) 



Figure 12. Inverse dynamic torque for non-zero I.C. (Practical approach) 
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inverse dynamic torque, assuming the initial flexible coordinate values at t=t 0 are all zeros. 
Case 4 : Initial errors in general cases (Practical approach) 

Since the initial conditions, which were used for trajectory planning and the inverse 
dynamic calculation, are estimated from the measured initial conditions before t=to*, and the 
flexible coordinate values are assumed as zero, it is natural to have initial tracking errors for all 
system states and an imperfect inverse dynamic torque profile. 

However, Figure 13 clearly shows the asymptotically converging trend of tracking errors, and 
no overshoot or vibration after t=tf\ 

Case 5 : Experimental results of general case 

The control scheme of the Case 4 has been implemented to the real experimental flexible 
manipulator using same conditions. As the end position sensor is not available, the end position 
behavior can be predicted from the joint angle and the strain measurements. Though the 
tracking performance in Figure 14 is not as good as the simulation result, the experimental 
results, which show no overshoot or residual vibrations, are good enough for contact or bracing 
control. Non-perfect tracking is presumed due to the joint friction because the coulomb friction 
force deteriorates the inverse dynamic feedforward torque profile especially at the low speed. 
Despite the initial tracking errors and the friction, the end point stops at the desired position 
without any overshoot. It indicates the nonminimum-phase system zeros are canceled by the 
inverse dynamic torque almost perfectly. Thus, this control method for moving initial conditions 
can be useful to avoid the large impact for contact or bracing. 

8. Discussions 

This paper presented a time-domain inverse dynamic method that calculate the torque 
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Figure 13. Case 4: Asymptotic tracking for general initial tracking errors 
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profile and the desired flexible coordinates. It also gives a new interpretation of the inverse 
system of flexible manipulators. The inverse system has been divided into causal, and anticausal 
systems. And, the dynamic system state space can be divided rigid, causal, anticausal system 
subspaces. 

The inverse dynamic torque can be interpreted in two ways, first, as a feedforward torque 
cancels the nonminimum phase system poles and zeros. As a result, the total system transfer 
function becomes a unit transfer function that makes the output exactly equal to the input. 
Second, we can interpret the inverse dynamic method as a input shaping technique. It functions 
like a noncausal notch filter that has zeros at the system resonant frequency. Therefore, the 
generated input torque profile doesn’t have the system resonant frequency content of the input. 

Many feedback control schemes require full state feedback, and assigning zero values to 
the desired flexible mode trajectories has been acceptable for the tracking control of flexible 
manipulators. Since the time-domain inverse dynamic method provides flexible mode 
trajectories that match the desired tip trajectory dynamically, it can be used for a trajectory 
generation method for many advanced feedback control schemes. 

If we can measure the I.C. and accommodate the I.C. to the desired trajectories 
generation, we will obtain perfect tracking with no initial tracking errors eliminating the transient 
period like the zero I.C. case. However, it is very difficult to accommodate the anticausal part 
I.C. to the desir-ed trajectory generation because of the backward integration of the anticausal 
system. This paper’s trajectory generation method with the fixed structure acceleration profile 
will generally have initial tracking errors because of the anticausal part. However, the tracking 
controller shows the fast tracking convergence in the simulation and experiments. At the final 
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stop, no overshoot nor residual vibration was observed, so that this tracking method will be 
useful to avoid the large impact for contact or bracing. 
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ABSTRACT 

When a machine tool is mounted at the tip of a robotic manipu- 
lator, the manipulator becomes more flexible (the natural frequencies 
are lowered ). Moreover, for a given flexible manipulator , its compli- 
ance will be different depending on feedback gains, configurations, and 
direction of interest. In this paper, the compliance of a manipulator is 
derived analytically, and its magnitude is represented as a compliance 
ellipsoid. Then, using a two link flexible manipulator with an abrasive 
cut-off saw, the experimental investigation shows that the chattering 
varies with the saw cutting angle due to the different compliance. The 
main work is devoted to finding a desirable cutting angle which reduces 
the chattering. 

1. INTRODUCTION 

1.1 Motivation 

In real world applications, robot manipulators are mechan- 
ically very rigid by design. This rigidity is necessary for high 
positioning accuracy; however, it becomes difficult to perform 
operations when a rigid manipulator contacts a workpiece. On 
the other hand, flexible manipulators can provide passive com- 
pliance due to their link flexibility With this structural com- 
pliance, certain applications such as cutting a workpiece can 
be performed with pure position control. Thus, the compliance 
can provide a simple, inexpensive solution for certain applica- 
tions that otherwise could not be achieved with position control 
alone [1,2]. 

Figure 1 shows the block diagrams of the overall architec- 
ture of a cutting process with pure position control An abrasive 
cut-off saw is mounted at the tip of a manipulator Its link flex- 
ibility is represented by a spring constant {Ki). The position 
feedback signal is measured at each joint. Due to the flexibility 
of the link, it is possible to regulate the force applied to the 
workpiece by controlling the position of the end-effector rela- 
tive to the workpiece. However, if the stiffness of the link is 
high, any uncertainty in the position of the workpiece, or er- 


rors in the position servo of the manipulator will induce very 
large cutting forces. Eventually these uncontrolled large cutting 
forces will shorten the life of the grinding wheel and the ma- 
nipulator. Also, the high stiffness (Ki) causes a high frequency 
osciilationor unstable chattering due to reaction forces from the 
workpiece. This behavior can be easily explained by a root lo- 
cus with increasing Ki assuming that the position-controlled 
robot is a linear mass- damper -spring system. Therefore, the 
compliance of the manipulator becomes one of the important 
parameters, and more compliance is desirable in the cutting 
process with pure position control. 

In this paper, the compliance of a manipulator is derived 
analytically, and its magnitude is represented as a compliance 
ellipsoid. It is shown that the compliance will be different de- 
pending on feedback gains, link flexibility, configurations, and 
direction of interest. Then, using a two link flexible manipula- 
tor with an abrasive cut-off saw, the experimental investigation 
shows that the magnitude of the chattering varies with the saw 
cutting angle due to the different compliance The final results 
show a range of cutting angles with acceptable behavior for a 
point in the workspace with a near circular compliance ellipsoid 



Figure 1. Block Diagram of Cutting Process 
with a Position-Controlled Manipulator 
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1.2 Chatter 

In the utilization of metal-cutting machine tools, vibra- 
tions are often encountered. The contact between the tools and 
workpiece gives rise to excessive variations of the cutting force 
which endanger the life of the tool. These vibrations belong 
to the class of self-excited vibrations. The source of the self- 
exciting energy is in the cutting process. Furthermore, in many 
cases, the self-excited vibrations are mixed with forced vibra- 
tions excited by various sources such as continuous spinning of 
the tools [3]. In this paper, the self-excited and forced vibrations 
are referred to as chatter. 

Considerable knowledge about the influence of kinematical 
parameters on the chatter has been assembled. As yet, however, 
neither a complete theoretical description and analysis has been 
accomplished nor reliable ways found for eliminating chatter in 
grinding [3,4]. By experimental trial and error, general guide- 
lines have been established to reduce the tendency for chatter. 
Among these are the use of soft-grade wheels, frequent dressing 
of the wheel, changes in dressing techniques, reduction in ma- 
terial removal rate and more rigid support of the workpiece [8] . 
Even though many parameters influence the chattering, this pa- 
per examines mainly the relationship between the cutting angle 
and the compliance of the arm. 

3. DYNAMIC MODEL 

3.2 Cutting Process 

Exact modeling of a cutting process can be very compli- 
cated [3, 4, 5, 6]. For simplicity, this paper assumes that the cut- 
ting forces consist of the normal cutting force which is in the 
direction of the approach angle of the saw and the tangential 
cutting force which is perpendicular to the approach angle. The 
relationship between the normal cutting force ( F n ) and tangen- 
tial cutting force (F t ) can be assumed to be Coulomb friction 

PJ- i e. 

^ = It = 0.3 ~ 0.4 

* n 

Notice that F„ is larger than F t . Also, F t and F n are a function 
of the depth of the cut. 

3.2 Flexible Manipulator 

Modeling a multiple link flexible manipulator is a compli- 
cated procedure. The deflection of the arm is approximated as 


a finite series of separable functions which are the products of 
mode shape functions ip^(x) and time dependent generalized 
coordinates q/tj(t): 


«i(*. t) = (*)<?/ ,;(0 


where i represents the link number and j represents the mode 
number . The equation of the flexible arm motion can be de- 
rived from several techniques, but the Lagrange’s formulation is 
known for its simplicity and systematic approach [9,10]. Using 
Jacobians to compute the velocity of a point, the kinematic and 
potential energies are obtained by integrating the velocity and 
position of the point over the total system. These energies are 
used in Lagrange’s equations. Therefore, the equation of the 
motion is 


'Mrr 

Mr,}( 


Mfr 

M/f J V 

it) L° 
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where q r contains the generalized rigid joint coordinates and qj 
contains the generalized flexible mode coordinates. M v is the 
inertia matrix and Ki represents the link flexibilities N r and 
N/ include nonlinear terms such as the Coriolis and centrifugal 
force in each coordinate. G r is the gravity force. T represents 
joint torques and F represents an end point external force when 
the contact with environments occurs. Finally, J r and Jj are 
the partitions of the Jacobian matrix for a flexible arm. 

4. COMPLIANCE OF ARM 

The equation of the dynamic motion for a flexible arm is 
obtained in equation (1). Since the tip of the manipulator pro- 
ceeds very slowly and in a small range during cutting, the mo- 
tion can be assumed to be quasi-static and linear. Ther efore, 
the equation of the motion is simplified to the following form 
assuming the acceleration terms and the velocity terms are neg- 
ligible. 



Figure 2. Definition of Cutting Angle 
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As shown above in the above equation (2), the gravity force 
also contributes to a stiffness force. If a joint angle PD controller 
is applied to the flexible arm, the joint torque T will be 

T = —KpAq r ~ K v Aq r 


where K p and K u are the feedback gains. This input torque can 
be interpreted as a spring and a damper force. Again q r can be 
neglected in quasi-static motion. If we combine all these forces 
for the stiffness matrix, 
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This stiffness matrix is always invertible. Therefore, it can be 
rewritten as • 


/ Aq r \ _ T5* + Kp 0 

U<// l ° ^ 



(3) 


Since this matrix shows the relationships between the end-point 
external force and the joint coordinates, it is necessary to change 
the joint space to the Cartesian space. Using the Jacobian re- 
lationship which is 

J /i(2j) (4) 


Substituting equation (3) into equation (4) will give 


AX = [Jr 



0 



(5) 


Since compliance is defined as ’displacement per input force’, 
we may say 


C=[J r Jf ] 



0 

K l 



( 6 ) 


This matrix, C, is called the compliance matrix for a flexi- 
ble arm with a joint angle PD controller. This compliance ma- 
trix includes not only link flexibilities, but also feedback gains, 
stiffness and configurations of the arm. The link flexibilities are 
represented by Ki and the feedback gains by Kp. The config- 
urations of the arm are represented by J r and G r . However, 
the compliance matrix does not incorporate the mass proper- 
ties of the arm which also may influence the arm’s behavior. If 
equation (6) is expanded, it can be rewritten as 

C — J r + (7) 

This form of the compliance matrix shows a major difference 
between a flexible arm and a rigid arm. In a rigid arm case, 
Ki is assumed to be very large. Therefore, we may ignore the 
second term of equation (7) although it is the dominant term 
in a flexible arm case. 

Since the compliance is represented with a matrix for a 
multiple link manipulator, various input force directions cause 
different directions and sizes of displacements. This may be 
explained with a linear algebra concept. Prom equation (5), 
we see that the compliance is simply a linear transformation 
that maps the end-point force F in R* into a Cartesian space 
displacement in ft? The unit sphere in i? 3 defined by 

F t F < 1 

is a mapping into an ellipsoid in R z defined by 



(,) Unit Input Force Sphere 00 Compliance Ellipsoid 
Figure 3. Compliance Ellipsoid 


AX T {CC T )~ x AX < 1 

This ellipsoid has principal axes A^, A 3 e 3 , A 3 e 3 where e l is a 
unit vector and A, is an eigenvalue of ( CC T ). We call this the 
’compliance ellipsoid’ Therefore, a unit force F in the direction 
a induces a displacement in the direction i3 as shown in Figure 
3. 


5. A CASE STUDY 

A large experimental arm designated RALF (Robotic Arm, 
Large and Flexible) has been constructed and is under computer 
control. RALF has two degrees of freedom in the vertical plane 
The length of each link is about 10 feet. At the tip of RALF, 
an abrasive cut-off saw is mounted as shown in Figure 4 Using 
the compliance ellipsoid, we explore the desirable cutting angles 
for RALF with acceptable chattering behavior. 

5.1 Analysis 

Based on the actual dynamic parameters of RALF, the dy- 
namic equation is derived in the form of equation (1). Then, 
actuator dynamics are assumed to be constant gains since their 
bandwidth is very high compared to the arm dynamics. The 
actuator gains are included in the feedback gains Kp. The 
nominal configuration during cutting is the following: the first 
joint angle is 106.6 degrees, and the second joint angle is 101. S 
degrees. The compliance matrix is computed and its magnitude 
is represented in R? with an ellipsoid. Simulation results show 
that the manipulator’s axis of least compliance is at an angle 
30 degrees with the horizontal and the axis of greatest compli- 
ance is at an angle of 120 degrees in Figure 5(a). Therefore, 
the 120 degrees cutting angle is desirable to produce the least 
chattering due to its greater compliance. Different shapes of the 
compliance ellipsoid can be obtained at different configurations. 
For example, if the first joint is at 110 degrees and the second 
joint is at 50 degrees, the compliance ellipsoid can be shown as 
in Figure 5(b). 

5.2 Experiments 

To measure the chattering in plane motion, two accelerom- 
eters are mounted at the tip of RALF . One accelerometer mea- 
sures the X direction vibration and the other accelerometer 
measures the Y direction vibration referenced to the manipu- 
lator base coordinates. Experiments use three different cutting 
angles 0, 40 and 90 degrees. The cutting angles are shown as in 
Figure 6 and the manipulator follows the given trajectory Each 
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Figure 4. RALF with an Abrasive Cut-off Saw 

case assumes that it has the same cutting parameters except for 
the cutting angle. For instance, cutting velocities are the same 
for each case, and the same feedback gains are used too. The 
workpiece is a half inch diameter steel bar and is much stiffer 
than the manipulator system itself. 

First, the abrasive cut-off saw moves very close to the work- 
piece. Then, the saw is turned on without contact with the 
workpiece. The vibrational signal is measured by a signal an- 
alyzer, and its power spectrum is plotted in Figure 7(a) . The 
first natural frequency is observed at 4.5 Hz compared to 5 Hz 
from the mathematical model. Also, another peak is observed 
at about 62 Hz This frequency is believed to originate from 
dynamic imbalance ofthe saw motor turning at 3800 rpm (63.3 
Hz) by the manufacturer’s data. 

Second, the cut-off saw followed the 0 degree desired tra- 
jectory by a joint angle PD control. The trajectory is computed 
based on Dickerson and Costing’s work [llj. It takes about 10 
sec for the saw to go through the workpiece. The acceleration 
power spectra are measured for 2 sec four times during cut- 
ting and are averaged to eliminate the influence of non-periodic 
noise. The same procedures are used for each experiment. In 
Figure 7(b), first peak is measured at 9 Hz. We may interpret 
this shift of the first frequency due to the change in the bound- 



Figure 5. Compliance Ellipsoid for RALF 
(a) when B x = 106.8 deg and 9? - 101.8 deg 
(b) when 0i =110 deg and 9 2 = 50 deg 


ary condition when the saw touches the workpiece. Also, we 
may notice that the rotation speed of the wheel is reduced due 
to the contact friction force. 

Third, when the saw cuts the workpiece at 40 degrees, the 
first natural frequency(9 Hz) no longer dominates as before and 
is mixed with other frequency signals in Figure 7(c). However , 
the higher mode at 37 H 2 becomes more noticeable. In other 
words, chattering becomes faster. 

Fourth, the 90 degree cutting shows smaller magnitudes of 
vibration in the Y direction , and the mode at 22 Hz becomes 
important (See Figure 7(d)). We expect that this angle will give 
the least chattering based on analytical analysis. Experimental 



workpiece 

Figure 6. Various Cutting Angle 
(a) 0 deg (b) 40 deg (c) 90 deg 
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Figure 7. Power Spectrum Measurement in X-Y direction 

(a) Without Any Contact with Workpiece 

(b) 0 degrees Cutting Angle 


Figure 7. Power Spectrum Measurement in X-Y direction 

(c) 40 degrees Cutting Angle 

(d) 90 degrees Cutting Angle 


data feils to show a distinct advantage. 

Finally, various cuts have been performed by a tele-operated' 
joystick under human control. Most of the cutting processes are 
successfully accomplished without any severe chattering. How- 
ever, its measurement is not included in this paper due to space. 
Experimentally, the contact velocity is one of critical factors 
which initiates chattering. 

6. CONCLUSIONS 

Taking advantage of the passive compliance of the flexible 
manipulator, certain applications such as cutting a workpiece 
are performed with pure position control. This provides a sim- 
ple, inexpensive solution for certain applications that otherwise 
could not be achieved with position control alone. 

Both computer-controlled cutting and human-operated cut- 
ting were performed with minor chattering. However, contact 
velocity should remain very small to reduce chattering. 

The contact with the workpiece causes a shift of the first 
natural frequency of a flexible arm due to the change of the 
boundary conditions. Different cutting angles produce different 
frequencies of vibrations due to the different compliances in the 


direction of forcing. 

Analytical studies predict 120 degrees cutting as the most 
desirable. However, this experimental investigation could not 
show distinct differences in the magnitude of chattering, al- 
though we may say that the 90 degree cutting angle generates 
a smaller magnitude of chattering in the Y direction. The com- 
pliance ellipsoid for our test bed is not elongated enough to 
make distinct differences in compliance. A different configura- 
tion could have made a more elongated ellipsoid, but further 
experimental tests have not been conducted due to physically 
limited location of the workpiece. 
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Abstract 

How do multibody systems move in free-fall? For instance, when 
a cat falls, it flips over before it reaches the ground. How does it do 
that? Multibody systems in free-fall move very differently than robots 
which are bolted to the ground. Once a robot with a fixed base stops 
moving, the link positions can be determined by kinematics alone. 
This is not true for a robot or multibody system in free-fall. The 
filial link positions of a robot in space depend on the link trajectories 
during the motion as well as its kinematics. Kinematics and dynamics 
are tightly coupled for multibody systems in free-fall. Given these 
difficulties, how can we plan motions for multibody systems in free- 
fall? 

The proposed research will center on several issues necessary to 
plan and execute motions for multibody systems in free-fall. 

L What motions are possible for a multibody system in free-fall? 
Mathematical techniques from nonlinear control theory will be 
used to study the nature of the system dynamics and its possible 
motions. 

2. How can we plan the link motions and joint torques necessary 
to move from one configuration to another ? Optimization tech- 
niques will be applied to plan motions. 

3. How can we store precomputed motion plans efficiently ? Since it 
is unlikely that motion plans can be computed in real time, pre- 
computation will be necessary. Image compression techniques 
are proposed to compress the precomputed motion data for stor- 
age. 

4. Once a motion is planned, how can the syste m execute the motion 
faithfully? A linearized controller will be devised to control the 
system while it executes preplanned trajectories. 

Symbolic manipulation techniques will be used in the research (where 
practical) to reduce chances for algebraic errors and to make the ap- 
proach easier to apply to new multibody systems in free-fall. 

The proposed research applies to a number of activities. Most 
obviously, it can be used to plan motions for robots in space. It can 
be used to plan limb motions to reorient astronauts. The research may 
also be useful to plan the movements of airborne divers, gymnasts, and 
jumpers. 
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Introduction 


1.1 Purpose 

The mot ion of ground-based robots is reasonably well understood. The un- 
derlying kinematic and dynamic analysis relies on the fact that robots are 
bolted to the floor — which will not move appreciably despite the motions of 
the robots. This is not true for robots in space. More generally, this is not 
true for multibody systems in free-fall. When a multibody system is not 
attached to the earth, its motions are considerably more complex than its 
ground-based counterparts. The dynamics and kinematics become' inextrica- 
bly coupled. 

Planning motions for multibody systems in free-fall is more difficult than for 
fixed-base robots largely because of the interaction of the kinematics and 
dynamics. Planning feasible or optimal trajectories will require extensive 
off-line computation and cannot be done in real time. A way is needed 
to precompute and store the possible motions so that the possibilities can 
be retrieved later and used quickly for real-time planning and execution 
purposes. 

The purpose of this research is to develop an end-to-end system that can be 
applied to a multibody system in free-fall to analyze its possible motions, 
save those motions in a database, and design a controller that can execute 
those motions. A goal is for the process to be highly automated and involve 
little human intervention. Ideally, the output of the system would be data 
and algorithms that could be put in ROM to control the multibody system 
in free-fall. 

This research applies to more than just robots in space. It applies to any 
multibody system in free-fall. This includes astronauts in space, falling mech- 
anisms, athletes in jumps or dives, and airborne gymnasts. 

1.2 Motivation 

To illustrate the complexities and potential of the types of motion that will 
be addressed by this research, consider the following examples. 
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Falling cat problem 


When a cat is held upside down and dropped, it manages to turn ilself 
rightside up before it lands. This well-known phenomenon has long intrigued 
children and scientists. Thomas Kane resolved the question in 1969 with a 
dynamic simulation that showed the motions necessary to turn the cat over 
during its fall [35]. His dynamic simulation involved only two bodies but 
duplicated the complex motion of the cat during the maneuver. 

Astronaut attitude change maneuvers 

Astronauts are taught a series of maneuvers that allow them to change their 
orientation while in free-fall. Again, it was Kane who developed these ma- 
neuvers [36. 37]. For each of the desired rotations (pitch, roll, and yaw), 
he developed simplified equations of motion and then analyzed them to de- 
termine what cycles of limb motions were necessary to give the desired ori- 
entation change. The result was several cycles of simple limb motions that 
produce orientation changes about the desired axis. 

Divers, Gymnasts, and Jumpers 

Some of the most complex motions of systems in free-fall occur when spring- 
board divers are in the air [21, 132, 133. 134. 135]. Gymnasts and ath- 
letes perform similar maneuvers while in the air during their activities. The 
movements of high jumpers while off the ground are also complex. The "Fos- 
bury Flop” revolutionized high jumping by improved maneuvers while in the 
air [125]. 

Robot servicing in space 

Using robots to service satellites in orbit is a goal of NASA [56. 101]. Several 
robot designs have been proposed which involve complex arms attached to 
large bodies with control-moment gyros and thrusters for maneuvering and 
station-keeping [S. 33. 57. 58], These approaches to motion control have their 
drawbacks. Control-moment gyros are complex and expensive. Thrusters 


produce plumes which may impinge on delicate equipment, 
could he developed to plan and execute attitude and configt 
via limb motions, simpler, safer, and less expensive servicing 
be designed. 


It techniques 
ration changes 
systems could 



2 Current Approaches (Literature Review) 


Many researchers have addressed various aspects of this research. Relevant 
research is reviewed below. 


2.1 Motion Planning and Control for Robots in free- 
fall 

Only a few researchers have directly addressed the problem of how to plan 
and control motions of a robot or multibody system in lree-tall. 

Kane’s work in astronaut maneuvering has already been mentioned. This 
research was reported in 1970-72 [3b. 37]. He worked out the equations ol 
motion for the human body and simplified them tor the desired rotations ol 
the body trunk. The resulting cyclic motions are interesting and useful but 
very different from the types of motion desired in this research. In any case, 
his work was specific to the human body and did not address the question 
of general configuration change of multibody systems in free-fall. 

Longman. Lindberg. and Zedcl considered a robot arm mounted to a satellite 
and developed special kinematics that dealt with the dynamics-kinematics 
interaction problem [57. 54. 58]. They assume the satellite base body con- 
tains reaction wheels to keep it from rotating when the arm moves. Their 
kinematics compensate for the translational movement of the base during 
movements of the arm. Their kinematic simplifications depend on the ab- 
sence of base rotation and do not generalize to multibody systems without 
reaction wheels. 

Vafaand Dubowsky developed the virtual manipulator technique for analyz- 
ing the kinematics and dynamics of robots in space [117. 118. 119. 120. 121. 
SO]. They devised a way to construct an imaginary manipulator with dimen- 
sions and inertial characteristics related to the actual system. The motion 
of the imaginary system and the real system are closely related and exactly 
the same for one point of the actual manipulator (such as the end efTector). 
Once the position of this common point is determined, the necessary virtual 
manipulator configuration is easily computed and then the corresponding 
joint positions of the actual system can be determined. The dynamics of 
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virtual manipulators also have the advantage that the conservation ot lin- 
ear momentum is implicitly integrated and eliminated from the equations 
of motion. Unfortunately, this approach does not provide any solutions for 
how to move from one configuration to another. In other words, the virtual 
manipulator approach can be used to determine the final joint positions to 
accomplish some task but is not very helpful in determining the necessary 
joint motions to move from the starting configuration to the final configura- 
tion. They did apply this technique to manipulator motion planning by using 
small cyclic motions to produce small motions of the manipulator base. Al- 
though this approach is useful for planning the motion of space manipulator 
end effectors, it has limited usefulness for planning large motions in which 
the entire final configuration is specified (such as in gymnastics). 

Given an end effector position that can be reached, the virtual manipulator 
approach can be used to determine the necessary joint angles of the space 
manipulator to acheive that position. 

Umetani and Yoshida studied continuous path control of manipulators. They 
devised an extended Jacobian for kinematic analysis of the motion of the 
end-effector [116]. The generalized Jacobian enforces the conservation of 
linear and angular momentum for the space manipulator. They use the 
new Jacobian to plan continuous motions of the end-effector in space and 
simulate them for a OMV (Orbital Maneuvering Vehicle.) This approach 
has similar benefits and limitations as Yafa's work because it concentrates 
on end-effector motion. 

Nakamura and Mukherjee also addressed the problem of planning motions for 
space robots [71]. Their work deals specifically with the nonholonomic nature 
of the conservation of angular momentum. They devise kinematics which 
incorporate the linear and angular momentum and of the space manipulator. 
They then devise a controller based on a Lyapunov function. Their approach 
is promising but initial results were disappointing because the controller could 
get stuck during the motion. They solved this problem in their later work 
[72. 73] by designing a bi-directional control algorithm. Unfortunately, the 
resulting motions involve unusual cyclic mot ions and other peculiarities which 
indicate that the motion is not very general. 

Sreenath and Krishnaprasad (among others) use mathematical approaches 
to attack the problem of controlling the motion of multibody systems in 
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free-fall [102, 103. 45. 131]. While these methods are diflinilt to understand, 
they do appear to offer promising techniques. Unfortunately, their research 
is not very useful to the proposed research because they currently consider 
only planar systems. In "Nonlinear Control of Multibody Systems in Shape 
Space," X. Sreenath indicates that extending their results to three space is 
“definitely non-trivial." [102. p. 1780]. 

In a recent paper[67]. Murray and Sastry use Chow’s theory and Lie brackets 
to determine whether a motion is possible for a system subject to a nonholo- 
nomic constraint linear in velocities. If the motion is possible, they construct 
a path using sinusoidal path segments. They apply their technique to the mo- 
tion of planar systems in free fall and to the nonholonomic vehicle problem. 
Their approach has useful insights and may be applicable to the proposed 
research but has not yet been applied to non-planar motion. 


2.2 Path Planning for Mobile Vehicles 

An area related to the current research is path planning for mobile robots. 
Much research has been devoted to planning land vehicle motion for indoor 
and outdoor vehicles. Works which considers mobile vehicles as points are 
not considered here since they are irrelevant to the proposed research. There 
are a few researchers who address the nonholonomic nature of vehicles with 
limited steering capabilities. 

Laumond considers a nonholonomic vehicle and proves that it is possible to 
plan collision free paths through a cluttered area by combining sets of small 
cyclic motions [46]. In later work he showed that whenever it is possible 
to plan a jagged path, it is also possible to plan a smooth path for the 
same motion [47]. Barraquand and Latombe addressed similar issues and 
devised planning techniques based on potential field techniques [5]. They 
apply their approach to difficult problems such as parallel parking a vehicle 
with several trailers. They also applied their approach to robot arms with 
many degrees of freedom. Similar research is covered by Jacobs and Canny 
[30, 31]. This type of research has many insights to offer for nonholonomic 
systems. Unfortunately, the nonholonomic constraints due to limited steering 
angles are simpler than the nonholonomic constraints due to the conservation 
of angular momentum. These approaches have not been applied to multibody 



systems in free-fall and it is not dear how applicable they are. 


2.3 Trajectory Planning and Control for Fixed-Base 
Robots 

An extensive amount of research lias been done on planning motions tor 
fixed-base robots. A sampling of this research is given in the references 
[6, 17. 34. 40. 22. 23, 7S. 64. 74. 75. 91, $9. 97. 98. 100. 105. 11 1. 126. 127]. 

One of the most promising approaches is presented by Tan and Potts [106. 
107. 108]. Their approach does everything. Their technique is intended 
for fixed-base robots but is general enough to handle multibody systems 
in free-fall. Their technique handles full dynamic nonlinearities, actuator 
limitations, joint constraints (position, velocity, and jerk), avoids obstacles, 
and incorporates an energy objective as well. This approach has not been 
adapted to multibody systems in free-fall but appears useful for this research. 

Another promising approach to planning opt imal motions is given by Luus in 
recent research on controlling chemical processes [61]. His approach is based 
on dynamic programming and may be useful for the problem at hand. Lmis 
has applied his technique to problems with up to eight nonlinear ordinary 
differential equations and determined optimal controls with limits on input 
variables. 


2.4 Tabular Planning and Control 

Another area relevant to this research is precomputing motion data and stor- 
ing for later use in planning and control. This is sometimes called a tabular 
approach since motion data are precomputed and stored in tables for later 
retrieval in planning or control. Very little has been done in this area. 

Raibert does use tabular techniques with some success for control of the 
cyclic parts of motion of his one-legged hopping machine [81, 82. S3]. Tab- 
ular techniques were also proposed by Albus [1, 2]. Hollerbach criticizes 
tabular approaches in general (and these in particular) when he concludes 
that dynamics simualtion codes can be made fast enough to run in real time 
[27]. This criticism is not relevant to the proposed research. It may be pos- 
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sible lo simulate the motion of multibody systems in live-lall faster that real 
time if the torques or forces to apply at each actuator during the motion are 
known before hand. The proposed research is to determine these actuator 
inputs. This cannot be done in real-time using known techniques even on 
super computers. 


2.5 Symbolic Manipulation 

Symbolic manipulation offers researchers many opportunities to improve 1 tin' 
quaility of their work by producing results much faster than is possible by 
hand, reducing the chance of mathematical errors, and allowing handling ol 
more difficult problems. Applying symbolic manipulation to robot kinematics 
and dynamics is not new. 

Hussain and Noble used symbolic computation for forward and inverse kine- 
matic analysis of specific robot geometries which assisted the user but still 
required considerable interaction[28]. Direls developed a system for manip- 
ulation of matrices with symbolic entries and used this to analyze robot 
kinematics [16]. Kircanski and Vukobratovic constructed a system using 
FORTRAN-77 to symbolically generate the forward kinematics and Jacobian 
of a robot but not the inverse kinematics solution [41]. Lloyd and Hayward 
applied MACSYMA to the same problem and derive kinematics for several 
common robot architectures [55]. Tunstel and Yira also use MACS5 MA to 
construct robot kinematics symbolically as an educational and design aid 
[113]. They also introduce a number of rules (symbolically implemented) 
that simplify the results. 

Many researchers have also developed dynamic equations of motion for multi- 
body systems symbolically. Liegois and company developed PL/1 software 
to derive equations of motion using a Lagrangian formulation. Others have 
written FORTRAN programs for symbolic generation of equations of mo- 
tion for multibody systems using various approaches: Newton-Euler [43. 44] 
and Kane's equations [IS]. Other similar work has been done by various 
researchers [7, 11. 24, 26. 29, 49. 50. 62, 66. 76, 77. 69. SS, 90. 94. 95, 96, 109, 
110, 111, 112, 122. 128. 129. 136. 137] Others have applied similar techniques 
to systems with flexible components [12. 59. 115]. Many of these systems gen- 
erate the equations of motion encoded in a FORTRAN or C program suited 
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to compiling and running for simulation purposes. Symbolic manipulation 
has also been applied to control applications [91. 10-1. 107]. 


2.6 Multibody Dynamics 

Multibody Dynamics is a huge field. Many people have developed widely 
varying approaches to the problem of simulating and controlling multibody 
systems. Several references cover Multibody dynamics in detail [-1. 19. 93. 
So. 130]. Others, too numerous lo mention, deal with dynamics in gen- 
eral and are applicable to multibody dynamics. Although serious multibody 
dynamics research was done more than SO years ago [20]. the field is noi ex- 
hausted. Recent developments include many recursive techniques lor inverse 
dynamics with operations counts proportional to the number of elements 
[3, 15. 19, 25. 27. 38, 39. 53. 60. 65. 86, 87. 124. 123]. (Most of these are based 
on recursive Newton-Euler approaches; some are based on operation space 
approaches [19. 53. 86. 87].) The most efficient of these approaches is given 
by He and Goldenberg [25]. Their recursive technique requires 9 1 ( // — 1 ) — 6 
multiplications and 86(n — 1) — 10 additions, where n is the number of bod- 
ies. The efficiency of these recursive techniques allows the computation of 
joint torques necessary to produce desired motions in real-time for reasonably 
complex svstems. Forward dynamics algorithms are not quite as efficient vet 
[48]. 
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3 Proposed Research 


Before analyzing the proposed research in detail, an overview may be helpful 
to orient the reader. 


3.1 Overview 

The goal of the research is to develop and test a system which can precom- 
pute, save, and execute motions for multibody systems in free-lall. 7 he basic 
components of the research are listed below. 

1. Analyze motion possibilities 

2. Implement simulation system 

3. Implement symbolic construction of equations of motion 

4. Design optimal controls to accomplish motions 

5. Implement symbolic generation of optimal control scheme 

6. Precompute motions between selected configurations 

7. Adapt compression techniques to compress motion data 
S. Design linearized motion tracking cont rol scheme 

9. Implement symbolic generation of linearized controller 

10. Use simulation to verify linearized controller 

11. AppI y system to example multibody systems 


3.2 Proposed Approach 
3.2.1 Definitions 

Several terms are used in specific ways in this proposal and are defined here. 
The terms appear in italics in the following definitions. 

Configuration (or post ) refers to the shape of the body as determined by the 
joint positions. Orientation refers to the attitude ol the system with respect 
to some global reference frame. More precisely, orientation refers to the 
attitude of some reference link of the body with respect to a global reference 
frame. If a robot is bolted to the floor, there is no reason to make the 
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distinction between configuration and orientation. Once the base of a robot 
or multibody system is free to move with respect to the global reference 
frame, this distinction becomes useful and important. 

Typical Configurations are configurations of the multibody system that occur 
often during motions of the system and are useful in studying and planning 
its motions. For instance, a tuck is a typical configuration for divers. For 
more detail, see Appendix A. page T2. 

Motions refer to movement from one combination of configurat ion and ori- 
entation to another combination of configurat ion and orientation. In this 
research, this will be accomplished strictly by joint motions. 

3.2.2 Analyze motion possibilities 

What motions are possible for multibody systems in free-fall? That question 
is central to this research. The possible motions depend on the nature of the 
mechanism, the initial configuration and orientation and the final configu- 
ration and orientation. For instance, a mechanism with one revolute joint 
(like a hinge) has a limited range of motion. It can open and close but the 
axis of the hinge cannot be tilted by opening and closing the hinge. This is 
because its motion is holonomic. A nonholonomic system has more potential 
motions. Consider a vehicle on the plane with a limited steering angle. The 
front wheel imposes a motion constraint that is nonholonomic. The vehicle 
has three degrees of freedom in the large but only two degrees of freedom at 
any instant. Yet, by careful maneuvering, any position in the plane can be 
reached. Multibody systems in free-fall must conserve angular momentum 
because they have no external torques acting on them. The conservation 
of angular momentum can be thought of as a nonholonomic constraint on 
the motion of the system in free-fall. Depending on the character of the 
angular momentum, a mechanism in free-fall may be able to move from any 
combination of configuration and orientation to any other combination of 
configuration and orientation: or it may not — as in the case of a hinge in 
free-fall. Obviously, since no external forces are used, the system center of 
mass will not move either case. 

This research will investigate this issue further and devise tests to be applied 
to determine if each of the desired motions is possible. For example. Frobe- 


ii ins' theorem can be applied using Lie brackets to evaluate the nonholonomic 
nature of the angular momentum (whether it is integrable) [51. 5. 71. 99. 70]. 
This can be done symbolically^ 2] since the angular momentum can be gen- 
erated symbolically. Research will also address the general controllability 
and reachability for these systems. It should be noted that it is very difficult 
to perform this type of research without symbolic manipulation due to the 
complexity of the equations of motion and angular momentum. 

3.2.3 Implement simulation system 

A basic part of the proposed research is a simulation environment in which 
the various components of the research will be implemented and tested. This 
simulation system will allow the user to const met robots from links and joints 
and then simulate kinematics and dynamics of the robots. The simulation 
environment will be used to verify the resulting motion libraries and control 
schemes. An extended description of how the simulation environment can be 
used is included in Appendix A. 

The initial implementation of the simulation environment will handle multi- 
body systems composed of rigid bodies since that is the focus of this research. 
To be even more useful, the simulation environment should also be able to 
handle flexible members. The software design and implementation will make 
provisions for future expansion in this direction. 

The software approach will be object-oriented and the code will lie written 
in an appropriate computer language such as C++, object-oriented Pascal, 
or Ada. An important component of such a system is the graphical display. 
These considerations and the goal of source-code portability indicate that 
C++ and X-Windows might be a good choice. 

3.2.4 Implement symbolic construction of equations of motion 

A number of systems exist for studying the motion of multibody systems. 
These include SD/Fast, SD/Exact. Autolev. DADS, and ADAMS. Others 
are mentioned in the literature review. These systems simulate multibody 
motions, and some generate C or FORTRAN code for simulation and con- 
trol purposes. Unfortunately, the output of most of these systems is not 
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directly suitable lor further symbolic manipulation (for controls analysis, for 
instance.) 

The proposed system will generate equations of motion in symbolic form 
suitable for further symbolic manipulation. (A few of the systems mentioned 
in the literature review do this.) The resulting symbolic form of the equations 
of motion will be used in three ways. First, the equations of motion will 
be used to generate executable code for simulation and planning purposes. 
Second, the equations of motion will be used to analyze system controllability. 
Third, the equations of motion will he used to const nut a linearized controller 
for trajectory tracking purposes. The last two will be done symbolically 
and the resulting symbolic material will be converted to appropriate code as 
necessary. 

The dynamical formulation that will be used has not been determined yet. 
A significant part of the research will involve comparing the various ap- 
proaches and choosing the most appropriate one to implement symbolically. 
Approaches to be compared include Newton-Euler. Lagrange equations (with 
Routh's extensions). Hamilton’s canonical equations. Kane's equations, and 
spatial algebra/screw theory approaches. 

To be suitable, the chosen technique of generating equations of motion must 
be suitable for symbolic implementation, and suitable for efficient simula- 
tions. The symbolic implementation should also apply typical methods to 
improve the efficiency of the code by doing such things as computing common 
subexpressions only once and by precomputing time-invariant terms. 

An issue to be addressed is how to adapt existing recursive approaches to 
multibody dynamics in free-fall. These formulations are satisfactory for sym- 
bolic manipulations for the systems under consideration. Unfortunately, the 
numerical implementations generally depend on the angular velocity of the 
base remaining zero. In free-fall, this is not true. In fact, the angular velocity 
and position of the base depends on the motion of all the joints (clue to the 
conservation of angular momentum.) Techniques to handle these unknown 
quantities during the recursions have not been described in the literature 
and will be studied in this research. This may lead to recursive formulations 
for the angular momentum of multibody systems in free-fall. (Note that 
this is not a problem for the symbolic use of recursive formulations since the 
unknown values are carried along as symbolic values in any case.) 
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3.2.5 Design optimal controls to accomplish motions 


Since the idea is to precompute motions, it makes sense to compute the best 
possible motions. Therefore, optimal controls approaches will be used to plan 
the motions. It should be noted that motions produced by this approach will 
be optimal in some sense but that the main goal is to plan motions that are 
feasible and avoid extensive cyclic motions. 

There are several possible approaches to be considered. An optimal control 
scheme based on variational analysis and the maximum principle is a logical 
candidate for computing the multibody motions. Appendix B gives a sample 
of the type of analysis proposed. The analysis actually used must satisfy 
several requirements. It must be simple enough and predictable to imple- 
ment via symbolic manipulation. It must produce the system of executable 
equations (for example, state and costate equations) which are reasonably 
efficient. The analysis in the appendix is given to illustrate the type of ap- 
proach proposed. 

Appendix C illustrates the application of the optimal control scheme from 
Appendix B to an example. 

Other approaches were mentioned in the literat ure survey. The technique de- 
scribed by Tan and Potts in "A Discrete Path/Trajectory Planner for Robotic 
Arms" is intended for fixed-base arms but is adaptable to this current prob- 
lem [107]. It involves constructing a discrete non-linear model of the robot 
dynamics[6S] which can then be used to construct a large non-linear pro- 
gramming problem. The approach is very flexible since it allows constraints 
on joint positions, velocities, jerks, and actuator limitations. It can avoid 
obstacles and will minimize a user specifiable cost function over the path. 

Luus has devised another technique based on dynamic programming which is 
also applicable. In the recent paper "Application of Dynamic Programming 
to High-Dimensional Non-Linear Optimal Control Problems." Luus used dy- 
namic programming to optimize several non-linear problems subject to input 
limitations. In one example, he studied a complex system with eight non- 
linear ordinary differential equations and determined optimum control input 
histories. 


3.2.6 Implement symbolic generation of optimal control scheme 

Once the optimal control scheme is designed, it must be implemented in 
terms of symbolic manipulations. An example ol the type implementation to 
accomplish this is given in Appendix C. In this example, the optimal controls 
scheme outlined in Appendix B is implemented. Sample results in terms of 
state and costate equations are given for example systems. F.xample code is 
also shown which has been generated from these state and costate equations. 

3.2.7 Precompute motions between selected configurations 

In order to prepare the system for movement between configurations in vari- 
ous orientations, the necessary motions must be precomputed. Tin 1 optimal 
control scheme must be applied to produce the movement data necessary lor 
each motion. This will be implemented in the simulation environment. 

The motion simulations will involve an extensive amount ol computation and 
may require assistance from fast mainframe computers. One advantage ol 
using X-Windows in a networked environment it is quite possible for the 
simulation environment to generate (' code for the motion simulation, move 
this to a remote computer (perhaps a super computer), compile the code on 
that computer, run it on that computer, and return the data to the simulation 
environment without user interaction. 


One premise of tins H’search is that the amount ol data generated by the 
motion simulations will not be over-whelming. To validate that assump- 
tion. it is necessary determine how much data will be generated for various 
situations. Appendix D contains a derivation of the number of data points 
that must be stored as a function of the various parameters. The resulting 
equation is: 

A dp = 3 A o A j A u A f, ( I ) 


where A op 
A o 
A 

A o 
A p 


Total number of data points required 
Number of relative rotations 
Number of joints or internal DOF 
Number of data points per variable 
(Number of time steps + 1 ) 

Number of poses (or configurations) 


To give some feel for the amount of data indicated by this equat ion, consider 
a few examples. Using two configurations and moderately optimistic values 
for the parameters in Equation 1. the amount of storage required tor several 
cases are given in Table 1 (see Appendix D for details). The first example 


Number of Joints . A j 

Aonr.-i (KB) 

3 

45.8 

6 

91.6 

14 

213.8 


Table 1: Amount of Data Necessary to Store Motions 

uses Nj = 3 and corresponds to a relatively simple mechanism. The second 
example uses A j = 6 and corresponds to a six degree- of- freedom robot. The 
last example uses A j = 14 and corresponds (roughly) to a human [134]. 

Although this is a large amount of data, it is within reason. Storing this 
amount of data on hard disks is quite feasible. Storing this amount of data 
in ROM is possible but A'p cannot be very large. 
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It may he possible to reduee the amount of data to lx* stored by only storing 
the joint positions during each motion. Joint velocities and joint torques can 
be computed on the fly by using recursive inverse dynamics formulations. 

3.2.8 Adapt compression techniques to compress motion data 

For each starting configuration and final configuration there will be three 
degrees of freedom in orientation that will be simulated. This can be thought 
of as a vector from the center of a sphere to some point on its surface plus 
an angle about that vector. One of the issues to be examined is how fine 
to subdivide these angles. The grid points must be close enough together 
so that interpolation between nearby motions will produce nearly correct 
results. This will be discussed further in the next section on the motion 
tracking control system. Unfortunately, increasing the number of divisions 
will tremendously affect the amount of number crunching necessary and the 
quantity of result ing data. 

Since the motion simulation will produce a. tremendous amount of data, an 
important component of this work will be how to compress it into a motion 
database (or library). Consider the plot for one joint position (or control 
input) over a motion. This is a single simple plot. Now consider a set of 
these for one of the degrees of freedom in orientation. Each plot of the joint 
position can be t reated like a scan line of an image so that the set of plots can 
be thought of as an image. There are three degrees of freedom so the resulting 
data can be thought of as a two dimensional array of images. Since the motion 
data can be thought of as images, one approach to compressing this data is 
to apply image compression techniques. The current state of the art in image 
compression for exact reproduction is roughly 3:1 for typical images. In this 
case, exact reproduction is not necessary: techniques exist which give nearly 
lossless single image compressions of roughly 10:1 to 40:1 [32. 10]. When 
a number of similar images are compressed, further compression is possible 
by exploiting the similarity between the images — resulting in compression 
ratios of up to roughly 100:1. With this type of compression, it is possible 
to compress an extensive set of data into a reasonable amount of space. It 
is expected that the set of pseudo-images will be relatively similar so that 
compression techniques will be effective. 
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The image compression techniques described typically depend on the im- 
age being composed of integer data with limited range, for instance. 0-2oo. 
The joint position and control input data will typically be floating point. 
An issue to be addressed is what level of quantization will allow acceptable 
reconstruction of the joint and control profiles. 

There may be a relationship between the type of compression scheme imple- 
mented and how the system will be used. If the mot ions are needed often and 
quickly (as it might be for planning), then retreival speed becomes an issue. 
The most effective image compression techniques depend on reconstructing 
t he entire image at once. All that will be needed in this case is the equivalent 
of one scan line from several different images. Some compression techniques 
may be more efficient for retreiving one scan line at a time horn an image 
(or set of images). 

The computations for compression will be extensive. This is not necessarily 
a significant problem for an actual application since video compression hard- 
ware exists today which can do the necessary compression at video frame 
rates. 

3.2.9 Design linearized motion tracking control scheme 

The process of compression means the reconstructed joint position and con- 
trol input profiles will not be exactly what they should be. Also, there will 
be uncertainties in the parameters of the actual system. Given the joint 
position and control input profiles necessary to accomplish some configura- 
tion and orientation change how can we persuade the system to actually 
complete the desired motion? Obviously some type of trajectory tracking 
controller will be necessary. There are a number of possibilities here. One 
is a time-varying linearized control system. Another approach is feedback 
linearization. In the research, various options will be examined. 

The controllability of the time-varying linearized system is an important is- 
sue that will be addressed. In a sense, the linearized control system controls 
the motions in the small at any instant. It will not be fully controllable (in 
general) since it cannot command motions that violate the angular momen- 
tum constraint. 
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A resonable approach (if the mechanism is suit al>l<-) is to reduce the order id 
the system used to determine the planned motion (for instance, by freezing 
some of tlx' joints). Then, during the motion tracking phase, the linearized 
controller can use those joints to keep the system close to the desired mot ion. 

3.2.10 Implement symbolic generation of linearized controller 

Once the form of the time- varying linearized controller is designed, it should 
not be difficult to use symbolic manipulation to apply it to the equations of 
motion. In this way. two implementation problems can lie addressed via sym- 
bolic manipulation. First, the system will generate C code to implement the 
linearized feedback control law. This code will not vary during the motion. 
Second, the system will generate C code to compute the time- varying data 
necessary for the linearized controller. This code will be run as necessary to 
update the data in the linearized feedback control code. 

3.2.11 Use simulation to verify linearized controller 

To test the motion data libraries and linearized controller, the system will 
use perform simulations. The system will use standard multibody simulation 
techniques with joint actuator inputs from the linearized controller and mo- 
tion libraries. These simulations will test many phases of the research. They 
will also give a feel for what kind of accuracy and resolution is necessary in 
the motion database to give adequate control with the linearized controller. 

3.2.12 Apply system to example multibody systems 

To illustrate use of the system and to test it. it will be applied to several 
example multibody systems in free-fall. Useful examples include two body 
systems, typical space robots, and simplified human models. Although hu- 
man motion in free-fall is a desirable application, it may be too ambitious 
for initial applications due to its large number of degrees of freedom. 
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3-3 Expected Results and Contributions 


It should be noted that no single piece of this research is revolutionary. 
At most small extensions from the state of the art are proposed. Y\ hat 
makes this research unique is the way the components are put together. 
Nobody has yet successfully addressed the end-to-end problem of how to 
control multibody systems in free-fall in real time. This will be the primary 
contribution of this research. 

Other contributions will include: 


• Extending recursive multibody formulations for numerical simulations 
of systems in free-fall. 

• Embedding the multibody tree structure in the soft ware objects created X 
to model it. 

• Construction of a flexible, powerful, and portable simulation environ- 
ment which can be applied to real motion problems. 

• Design and implementation of optimal control for configuration change. 
Using symbolic manipulation to construct the optimal controller. 

• Using image compression techniques to compress motion data. 

• Storing precomputed motions for complex systems for later use. 7 

• Using symbolic manipulation to implement the time- varying linearized 
motion tracking controller. 

• Use of symbolic manipulation for dynamics and controls in one inte- 
grated system. 
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Appendix A. How the Simulation Environ- 
ment Might Be Used 

The following description explains how the simulation environment might be 
used in a step-by-step manner. This secpience described here is not the only 
way the system ran be used, but does illustrate the basic components of the 
simulation system. 

1. User constructs system in simulation environment. The start- 
ing point of analyzing the motion of the multibody system is lot the 
user to model the multibody system to be analyzed in the simula- 
tion environment. This could be done by direct manipulation (on i he 
computer screen) or by reading an appropriate data file. To specify 
the multibody mechanism by direct manipulation, the user will select 
the links from a catalog of link shapes, specify (or modify) the link's 
geometric and inertia properties, indicate where any joints would be 
located, and what other links are attached to each joint. The system 
would then create software objects to model the links and joints. Note 
that this would automatically establish the connectivity (or tree struc- 
ture) of the system. The software for each object would know how to 
construct the relevant transformations (symbolically and numerically) 
to determine the robot kinematics. Similarly, the object's code would 
also know how to add their components (symbolically and numerically) 
to recursive dynamic formulations. 

2. User chooses typical configurations At this point, the user will 
choose typical configurations or poses tor the multibody system. The 
typical poses will be selected to put the multibody system in vari- 
ous useful or desirable configurations. For example, it the system is 
a human, a typical configuration might be a straight body with arms 
extended. In [21], Cliff Frohlich gives nine different human body config- 
urations that are commonly used by divers and trampolinists. Whether 
the body is upside down or right side up is not important in specifying 
the configuration. The typical poses will probably include only con- 
figuration information (such as joint positions) and will not include 
joint velocities. Including initial and final velocities in the typical con- 


figurations will increase computation and storage requirements by an 
unreasonable degree. The system will be able to plan individual mo* 
tions with initial and final velocities but these will not be included as 
part of the precomputation part of the research. This may limit the 
usefulness of the precomputation and storage aspect of the research to 
human motion since large initial and final velocities are often part of 
at li let ic mot ions. 

The typical configurations may also be selected to simplify the system 
dynamics. For instance, the pose might pul the wrists into a neutral 
position. During the planned motion, these joints might not be used 
to reduce the dimensionality of the problem. 

3. System automatically constructs equations of motion. Alter 
the description of the system is entered, the system will generate the 
equations of motion in symbolic form. The simulation system will be 
able to deal with the dynamics of the system in three ways. First, 
it will be able to simulate the multibody system dynamics directly 
using standard recursive approaches. Second, it will be able to generate 
the equations of motion in a symbolic form. Third, it will be able to 
simulate the multibody dynamics by using software code generated 
from the symbolic representations of the equations of motion. 

4. System automatically generates optimal controls. Once the 
equations of motion are generated in symbolic form, the optimal con- 
trol law for reconfiguration will be derived symbolically. This is why 
it is important to generate the equations of motion in symbolic form. 
This will also include generation of code to verify whether the motion is 
possible (from the analysis of nonlinear controllabiltiy and reachability 
analysis). 

5. System simulates optimal motions. Once the equations of motion 
have been generated and the optimal control scheme has been con- 
structed. these will be used to simulate optimal motions for orientation 
changes between the typical configurations. During each simulation 
the basic data of control inputs and states during the motion will be 
saved. The goal is to simulate the motions from any configuration to 
ariy other configuration in any other orientation where such motions 
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are possible. Hopefully, the number of typical configurations will be 
small. (If the number of typical configurations is large, the amount ol 
computation and resulting data may be excessive.) 

(i. System compresses the resulting simulation data into library 
of motion data. The preceding step will generate a large amount 
of data. The resulting motion data will compressed using image com- 
pression techniques to construct a database o( motion data (or motion 
library). 

7. System designs a linearized controller to execute the maneu- 
vers in the motion library. Alter the' system constructs the motion 
library, the library can be used to plan and execute motions. How- 
ever. the data compression techniques will introduce some errors into 
the reconstructed motion profiles due to quantization and other effects. 
This, along with imperfect modeling, indicate a motion tracking con- 
troller will be necessary. At this point, the simulation environment will 
use the symbolic version of the equations of motion to symbolically gen- 
erate the necessary linearized controller to allow the system to execute 
a ret reived motion profile. This linearized controller will be converted 
from symbolic form to usable software code for simulation purposes. 

S. User uses system to simulate various motions. Once the previous 
steps have been completed, the simulation environment can be used to 
simulate motions and test the motion tracking controller. This could 
be the goal of the entire system. Consider how such a system could be 
used: 

• The user could simulate possible motions just to see what they 
look like and what types of control inputs are necessary. 

• The system could be used to verify the linearized controller by do- 
ing simulations with precomputed motion data and a perturbed 
system. The control inputs could go into a true dynamic simula- 
tion to verify the results. 

• A diving coach could use the system to construct a new dive se- 
quence and show it to divers in a movie lorni. This would in- 
volve using several intermediate poses and splicing together the 
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necessary motions. (It might also involve general ing specilic new 
motions or new motion libraries.) The' system would take the sets 
of joint position profiles and construct a longer sequence to show 
how the motion would look. Individual joint motions could be iso- 
lated to show the diver what to do and when. It is even possible 
that this system could be used to discover maneuvers that have 
newer been thought of before. 

• Similar techniques could be used with astronauts to train them to 
do combined configuration and orientation maneuvers. 

• The system could be used to construct a motion library and track- 
ing controller for an orbital servicing robot. Ibis might involve' 
generating new motion libraries tailored to specific tasks. Tin' 
resulting data and code could be put into ROM for use' on orbit. 

One of the goals of this work is to reduce the amount of user interaction 
necessary. Ideally, steps 1 and 2 would be the only steps the user would have 
to supervise. In reality, some interaction will probably be necessary during 
some of the other steps as well. 


Appendix B. Sample Optimal Control 

Analysis 



Suppose that we have a nonlinear system such as a mullibody robot in free- 
fall and we would like to choose a set ot control inputs to move the system 
from one configuration to another. Start by putting the equations of motion 
of the system into tin* following form (this can be done symbolically Irom the 
equal ions of mot ion ): 

4x = f(x) + g(x)u 12) 

of 

where x is a vector describing t he st ate and velocities of 1 lie mult ibody system 
and u is a vector of available control inputs (such as control torques at each 
joint). This form was chosen since the equations of mol ion lor multibody 
systems can generally be put into this form. 

The initial configuration. x(/ 0 ) = x 0 , is known and the goal is to use the 
control inputs to move the system into the final configuration, x(//) = X/. 
(where tj is unspecified). The requirement that the system acheive the desired 
terminal state can be formulated in the terminal constraint: 

\P[x(t/). //] = x — X/ = 0 (3) 

The terminal constraint is adjoined to the terminal cost (which is zero so tar) 
by using a vector of constant multipliers. v\ 

4>[x(//). tj] = v 7 # = i/ r (x - x/) (4) 


The motion should minimize some combination of time required for the mo- 
tion and control effort of the actuators during the motion. A suitable cost 


function is: 


where: 


: $[X( </M/] + ! ’ L(ll 

Jt o 

(5) 

1 r 

L = a + -ir B u 

(6) 

0 < a < 1 

(T) 

f 6,(1 -«) > 0 if i = j 

U,J "10 if i* J 

(8) 
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The constant a is the ratio between the con Hiding goals of minimizing the 
time required for the motion and minimizing the control effort ol the actua- 
tors during the motion. If o = 0, then time is of no concern. If a = 1 then 
control effort levels are of little concern. Note that B is positive definite so 
Z? -1 exists. Also, since B is diagonal. B~ l = B~ T . 

Following the typical optimal controls approach, the previous constrained 
problem can be converted to an unconstrained optimization. I his is done 
]>v constructing a modified cost functional which enforces the equations ol 
motion by adjoining the equations of motion with lagrange multipliers. A: 

.7 = 4>[x( //)•//]+ f [A + A'f-x + f + gu)j <li OH 

To simplify the problem, the Hamiltonian of the system at some instant in 
time is introduced: 

H = L + A T (f + gu) ' (10) 

Which means the modified cost functional is: 

3 = + r [// - A r x] ,H (11) 


Taking the total variation of .7 and integrating by parts results in: 




Note that the variation with respect to A has been omitted since it leads 
back to the equations of motion. To force the variation of the modified cost 
functional. 7>J. to vanish, we choose: 


A 7 

= - 

~ (13) 

OX 

A T (//) 

= 

’(W 

ebe 

(14) 

OH 

c)u 

= 0 

do) 

ta| <= , 0 

= 0 

(16) 
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Those diokvs guarantee that the resulting A and u produce a stationary 
value of j. For J to be minimized, the second gradient must be positive 
definite: 


Hr* Hra 

IT U 

** u.r * * it u 


> 0 


\vher«* t he subscripts represent partial derivatives. This can be implemented 
symbolically. If this is satisfied, then the resulting controls will be optimal. 
See [9. p.oO]. for a detailed derivation of this requirement. Note that II „„ = B 
which was chosen to be positive definite. 


Applying these results to the problem at hand gi\(*s: 


Oil _ ■ T 

Ox ~ 




=> u T B + A r g = 0 


(IS) 


(19) 


Equation ( IS) gives the differential equations lor the costate vector. A. Solv- 
ing Equation (19) for u gives the optimal control law: 


u = -B~ l g T X 


(20) 


Since the terminal time, i j. is not specified, it is a free parameter. Treating / / 
as a free parameter produces a modified total variation of the cost functional. 
8J. The previous analysis and choices force all the terms to vanish except 
for the variation due to possible changes in the final time: 

btj ( 21 ) 

'='/ 


6J = 


0<t> „ 

w + w 


See [9. p.72]. for a detailed derivation of this requirement. Since $ is au- 
tonomous. = 0 and therefore Hj — 0. 

The Hamiltonian. H. is constant as can be seen by taking its derivative with 
respect to time: 



OH. OH- OH. OH 

— : — x + -rr-A + — u + — — 
Ox OX Ou Of 


( 22 ) 
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OH 

Ox 


f + f 


7 


OH 

Of 



+ 0u + 


OH 

~dT 


(23) 


(24) 


Tlie Hamiltonian is autonomous so jfjH = = 0. This means that H is 

constant on the optimal trajectory. Its constant value must he the same as 
its final value. Hj = 0. Therefore. If = 0 on the optimal trajectory. 

The final value of the costate vector. A. is determined from the' terminal cost : 


A(//) = A / = 


0<t> 

Tix m 


t=t f 


U T V)\ 


\ t=t f 




J t=tf 


= u 


(25) 

(26) 

(27) 

(28) 


This indicates that the final value of the costates are unknown constants. In 
order to determine their values, consider the Hamiltonian at the final time. 
Substitute u = -5 _l g 7 . A = u. f/ = f(x(?/)). and g / = g(x(//)) into H\ 


H = 


a + -u r Bu + A 7 (f + gu) 


(29) 


'(=</ 

-l„f. 


= o + -{-B 'gji/) 7 B(-B ’gji/) + ({f + gf(-B %jv)) (30) 

= a + + v T fj - (31) 


= o + i/ 7 f/- -i/Mg fB 'gj)v 


(32) 


At the final time. H is: 

H = a + u T ( ) - \u T {gjB~ x $)v = 0 


(33) 
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However, since II = 0 during the entire mol ion. ;i similar stcilement is t rut' 
at the intial t ime: 


H = o + fy — 


^H T {g 0 13 'g 0 r )Ai = 0 


(34) 


where n = A(/ 0 ). f 0 = f(x(/ 0 )). and g 0 = g(x(/ (J )). Not.- that n is not known 
hut that fu and gu are botli known. (Not.- tliai if ihis is satisfied at the 
initial time, the optimal control guarantees that the similar requirement will 
he satisfied at the final time.) 


This is a quadratic form in /i. It describes a multi-dimensional surface. 
Mathematically, this equation can have either no sob it ions, a unique solid ion. 
or many solut ions. 

A Newton-Raphson style scheme can be used to find a solution lor n if one 
exists. Consider the change in II to a small change in //, near a solution: 


H 0 (Hi.H2 /«. + A//, •.//,•+ 1 /'») 


Haiti) + ^± t ', + 0(tfm 

Ol>, 


If fi is near a solution, then Pi + A//,.//,+i /'„) ~ 0 and 

the 0{pf) terms are nearly zero so the resulting equation is: 

0 ss H 0 (n) + — (36) 

Op, 


This equation can be solved for A p,: 



(37) 


Then can be improved: 


(/h)neir — T A/// 


(38) 


This is a scalar equation for one //,. The values of : j~f- can be determined as 
follows: 

“ /* (go B 'g' ) (39) 

For computational efficiency, the n update equations of the form of Equa- 
tion (38) can be applied in parallel. 


Some experimentation may lx* required to choose llie weighting values a ami 
/>, to produce reasonable trajectories with reasonable joint actuator levels. 

In summary, the system has a set of first-order ordinary differential equations 
for state and costate. It is a two point boundary problem where initial and 
final values exist for the state. The initial values of the costate are unknown 
but candidates can be found using the Newton-Raphson type approach just 
described. This problem can be attacked by a shooting technique. Errors at 
the end can be used to improve the guess ol the inital costale values. The 
only moderation needed to this is that the initial costate values must be 
further modified so that they satisifv the quadratic form. 
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Appendix C. Applied Optimal Controls 

Example 

Appendix C.l. Description 

The purpose of this example is to demonstrate the type ol symbolic ma- 
nipulation techniques that ran be applied to generate optimal controls laws. 
Obviously, t his requires t hat the equal ions of mot ion are available in symbolic 
form. 


Appendix C.2. MACSYMA Usage Descriptions and 

Code 

The following description of the function OPTCONT is from the 
file OPTCONT. US AGE: 

The function OPCONT will take an array of first order ordinary 
differential equations with a cost functional and it will derive 
and return the optimal control, the costate equations, the 
Hamiltonian, the Hessian of the system. To use this function, the 
dynamic system must be put in the form of a list of first order 
differential equations of the form: 

dx 

— = f(x,u) 
dt 

USAGE: OPTCONT (odes ,L ,x_name ,u_name) ; 


INPUTS: 


odes : list of first order ordinary differential equations, eg, 
[dxl/dt = fl(x,u), dx2/dt = f2(x,u), . . .] 
describing the dynamics of the system. 

L : the cost functional of the system (in terms of x and u); 
[scalar function of x,u] 
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x_name : list of names of the states used in odes and L. 

Must be in the same order and in one-to-one corrspondence 
with the variables in the left hand sides in the odes. 

u_name : list of names of the control inputs used in f and L. 


OUTPUTS : 


A list composed of: 

costate : List of the costate first-order ode equations in terms 
of the original state variable names and new multiplier 
variables, Li. 

costate_names : a List of the newly introduced costate names 

u_opt : List of Optimal controls 

H : The Hamiltonian of the system, H = L + LT*f [scalar 
function] 

H.hess : The hessian of H, 

H_hess = [ Hxx Hxu ] 

[ Hux Huu ] 

This can be used to determine if the generated control 
inputs are optimal, [(n+ra) x (n+m) MATRIX of scalar 
functions] 

ROTE: The list costate contains variables of the form L1,L2,L3... 

This function uses KILL on all of these it uses, so existing 
variables with names of this form will be destroyed. 

By: Jonathan M. Cameron 

The MACSYMA code for the function OPTCONT follows: 

OPTCONT ( odes , L , x_name , u.name ) : = 

BL0CK([ n : Length(odes) , /* Humber of states ♦/ 

m : Length(u_name) , /* Humber of controls */ 

Lambda , U_Eqns , Hx , Hu , i , j , 
costate, u_opt, H, H_hess] , 
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/* Construct list of lambdas */ 

Lambda : □ , 
for i : 1 thru n do 

Lambda : append ( Lambda , [cone at ( 'L , i)] ) , 

Apply ( 'KILL, Lambda) , 

/* Form the Hamiltonian */ 
fl : L, 

for i : 1 thru n do 

H : H + LambdaCi] *RHS(0des[i] ) , 

/* Construct costate equations */ 

Hx : □ , 

for i : 1 thru n do 

Hx : append(Hx, [diff (H,x_name[i] )] ) , 
costate : [] , 
for i : 1 thru n do ( 

Depends (LambdaCi] ,T) , 

costate : appendC costate, [diff (lambda [i] ,T) = -Hx[i]]) 

), 


/* Solve for the optimal controls */ 

Hu : □, 

for i: 1 thru m do 

Hu : append (Hu, [diff (H ,u_name[i] )] ) , 

U.Eqns : □ , 
for i : 1 thru m do 

U_Eqns : append ( U_Eqns , [0 = Hu[i]]), 
u_opt : Solve (U_Eqns , u_name ) , 

/* Generate the Hessian */ /* H_hess = [Hxx Hxu] */ 

H.hess : ZEROHATRIX(n+m,n+m) , /* [Hux Huu] ♦/ 

/* Do the Hxx block */ 
for i : 1 thru n do 
for j : 1 thru n do 

H_hess [i , j] : diff(Hx[j], x.name[i]), 

/* Do the Hxu block */ 
for i: 1 thru n do 
for j : 1 thru m do 

H.hessCi, j+n] : diff(Hu[j], x_name[i]>, 

/* Do the Hux block */ 
for i : 1 thru m do 
for j : 1 thru n do 

H_hess [i+n, j] : diff(Hx[j], u_name[i]), 
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/* Do the Huu block */ 
for i : 1 thru m do 
for j : 1 thru m do 

H_hess [i+n, j+n] : diff(Hu[j], u_name[i]), 

/+ Return the results */ 

[costate, lambda, u_opt , H, H_hess] 

)$ 

The following description of the function SIMPCONT is from the 
file SIMPCONT. US AGE: 

SIMPCONT is a function which will take lists of state and costate 
first order ordinary differential equations, the optimal control, 
and other outputs of OPTCONT and will substitute in the optimal 
control in the state DDEs and then solve as many of the ODEs as 
possible. The resulting system is ready to simplify via boundary 
conditions . 

USAGE: 


SIMPCONT(state_odes , costate_odes ,u_opt , state.names ,costate_names ) ; 


INPUTS: 


state_odes : List of first order ordinary differential equations 
[dxl/dt = fl(x,u), dx2/dt = f2(x,u), . . .] describing 
the dynamics of the system. 

costate.odes : List of first order ordinary differential equations 
derived for optimal control by OPTCONT 

u_opt : list of optimal controls derived by OPTCONT 

state_names : list of the names of the states 

costate_names : list of the names of the costates (generated by 
OPTCONT) 


OUTPUTS: 



A list composed of: 


new_system : List of the new simplified set of state and costate 
equations 

new_names : List of the names of the states or costates in 
new_system 

constants : List of the new constants generated by solutions of 
ODEs 


By: Jonathan M. Cameron 

The MACSYMA code for the function SIMPCONT follows: 

SIMPCOKT( state , costate ,u_opt , state_name , costate.name) : = 

BLOCK ( [ n : LengthC state) , 
m : Length (u_opt) , 
new_system : □ , 

new_name : [] , 

const.num : 0, 
constants : □ , 

i, ii, soln], 

/* Substitute the optimal controls into the state equations */ 
for i : 1 thru m do 

for ii : 1 thru n do 

stateCii] : LHS(state [ii] ) = subst (u_opt , RHS(state [ii] ) ) , 

/* check each of the costate ODEs and try to solve them */ 
for i : 1 thru n do ( 

soln : ode2(costate[i] , costate.name [i] ,t) , 
if soln # 'FALSE then ( 

soln : subst (concat( 'C, const_num) ,7*C, soln) , 
constants : append ( [concat( 'C, const _num)] , constants), 
const_num : const_num + 1 , 

/♦ Do substitutions with soln to eliminate the costate */ 
for ii:l thru n do ( 

state [ii] : LHS(state [ii] ) = subst (soln, RHS(state [ii] )) , 
costate[ii] : LHS(costate [ii] ) = subst (soln, RHS(costate [ii] )) 
), 

for ii: 1 thru length(new_system) do 
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new_system[ii] : 

LHS(new_system[ii] ) = subst (soln, RHS(new_system[ii] ) ) 

) 

else ( 

new_system : append (new_system , [costate [i]] ) , 
new_name : append (new_name , [costate_name[i]] ) 

) 

), 

/* check each of the state ODEs and try to solve them */ 
for i: 1 thru n do ( 

soln : ode2(state[i] , state.name [i] ,t) , 
if soln # 'FALSE then ( 

soln : subst ( concat ( 'C, cons t.num) soln) , 

constants : append ( [concatC * C , const _num)] , constants), 
const_num : const_num + 1 , 

/* Do substitutions with soln to eliminate the state */ 
for ii: 1 thru n do ( 

state[ii] : LHS(state [ii] ) = subst (soln ,RHS(state [ii] ) ) , 
costate [ii] : LHS ( costate [ii] ) - subst(soln,RHS(costate[ii] ) ) 
), 

for ii : 1 thru length(new_system) do 
new_system[ii] : 

LHS(new_system[ii] ) = subst (soln, RHS(new_ system [ii] )) , 

/* Add this solution to the system */ 
new_system : append (new^sys tern , [soln]), 
new_name : append (new_name, [state_name [i] ] ) 

) 

else ( 

new_system : append (new_sys tern, [state[i]]), 
new_name : append (new .name, [st ate .name [i]] ) 

) 

), 


/* Return the results +/ 
Declare(constants, constant), 
[new_system, new_name, constants] 


)$ 
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Appendix C.3. Sample MACSYMA Session Output 

The following output is from a MACS') M A session using the functions OP T- 
CONT and SIMPCONT on a simple problem. 

(C3) load(optcont) ; 

Batching the file USERD1 : [CAMERON . PROP]optcont .mac ;60 
Batchload done. 

(D3) USERD1 : [CAMERON . PROP] optcont .mac ; 60 


(C4) load(simpcont) ; 

Batching the file USERD1 : [CAMERON . PROP] simpcont .mac; 33 
Batchload done. 

(D4) USERD1 : [CAMERON . PROP] simpcont .mac ; 33 


(C5) kill(x,v,u) ; 

(D5) D0NE 


(C6) dependsC [x,v,u] ,t) ; 

(06) DCCD, V(T), U(T)] 

(C7) state_eqns : [diff (x,t)=v,diff (v,t)=u] ; 

dX dV 

(D7) [" = V. - = U] 

dT dT 


(C8) state_names : ['x,’v]; 

(D8) £*« 
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(C9) control .names : [*u]; 

(D9) [U] 


(CIO) L : 0. 5*u~2; 
(DIO) 


2 

0.5 U 


(Cll) /* Find the optimal control and costate equations */ 
results : optcont (state. eqns ,L , stat e.names , control_names)$ 


(C12) costate.eqns : results Cl]; 

dLl dL2 

(D12) C = 0, = - LI] 

dT dT 

(C13) cost ate. names : results [2]; 

(D13) [LI, L2] 

(C14) opt. control : results [3]; 

(D14) [U = - L2] 


(C15) H : results [4]; 
(D15) 


2 

LI V + 0.5 U +L2U 


(C16) Hessian ; results [5]; 

[ 0 0 0 ] 

[ ] 

(D16) [ 0 0 0 ] 

[ ] 

[ 001 ] 


49 



(C17 ) /* Solve as many of the ODEs as possible */ 
results : 

simpcont (state.eqns , costate_eqns , opt.control , state.names , c os tate .names )$ 
DUBO: [MACSYMA.412. 0DE]ode2 .f as ; 1 being loaded. 


(C18) system : results [l] ; 

2 

CO 1 

(D18) 


CO T CO T 

[X = T ( Cl T + C3) + C2 , V = Cl T + C3] 

2 2 


(C19) names : results [2]; 

(D19) CX, V] 


(C20) constants : results [3] ; 

(D20) [C3 , C2, Cl, CO] 

(D21) DONE 

The next step depends on the problem to be solved. In this case, it is not hard 
to apply inital and final state values to resolve the resulting constants. In this 
example, the costate and state differential equations were solved completely. 
This will not happen with the type of systems to be considered in this re- 
search. In general, some differential equations will be produced. In any case, 
it is not difficult to take the results of SIMPCONT and use MACS'! MA 
to convert it to C or FORTRAN code {or simulation purposes. Symbolic 
manipulation systems such as MACSYMA and Mathematics have powerful 
capabilities to generate program code. A function to take the results of SIM- 
PCONT and generate code could also perform various optimizations such as 
computing common terms only once. 
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Appendix D. Movement Library Size 

Requirements 

In order lo validate the premise that the amount of data that will be saved 
is not too excessive, an estimate is presented in this appendix. 

The motions will move the system from one combination ol orientation and 
configuration to another combination ol orientation and configuration. Sup- 
pose the goal is to move from one typical configuration to another typical 
configuration. Since the system is in free-fall, the final configuration has 
three degrees ol attitude Ireedom with respect to tin* starting configmat ion. 
Think of this as the points on a unit sphere and another degree ol Ireedom 
about a line from the center of the unit sphere to the points on the surface 
of the sphere. In order to tesselate the unit sphere, a procedure based on 
constructing geodesics from icosahedrons can be used [13]. The degree ol 
the tesselation is Q. To tesselate the unit sphere so that there is an angle of 
approximately o between vertices, and: 

Q = int [arctan(2)/o] HO) 

where Q = degree of tesselation (41) 

n = approximate angle between verteces 

Therefore, the total number relative orientations that must be considered for 
moving from one configuration to another is: 


A'o = Number of relative orientations (42) 

= ^[lOQ 2 + 2] (43) 

For each orientation, the joint positions, velocities, torques must be saved 
over the motion. So the number of data points for one motion is: 

Xm = Number of data points per motion (44) 

= 3 XjX D H5) 

where Xj = Number of joints (46) 


\ D = Number of data points per variable 
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Multiplying ,\o and A \/ gives the mi tuber of data points necessary to store 
the motions from one configuration to another. 


A.\//> = Number of data points to move from (17) 

one configuration to another (-IS) 

= AoA,\/ (49) 

= 3A(;AjA/5 (70) 


Now assume there are A r typical configurations. The total number data 
points to for motions front any configuration to any other is: 


A op 


Total number of data points 



3AoAjA/j [Ap( Ap — 1 ) + Ap] 
3A'o AjApA p 


( 31 ) 


■>j 


(A3) 

(54) 


where the term (‘V) in Equation (52) gives the number of pairs of config- 
urations. The factor of 2 is necessary because the motions for each pair of 
configurations could considered in either direction. The next term. A p. ac- 
counts for motions from one configuration back to the same configuration in 
a different orientation. 


Given a data compression ratio of C’r to 1 and assuming that the floating 
point value for each data point can be quantized to 8 bits (or a byte), the 
amount of data (in KB or kilobytes) is: 


i DATA — Total amount of data in K Bytes 

A DP 


10240 


KB 


(55) 


(56) 


To give some feel for the amount of data indicated by these equations, con- 
sider a few examples. For 10° orientation resolution. Q ss 6 and A'o = 13.032. 
For each example A p = 10. Ap = 2. and Cr = 100. The results are given in 
Table 1 on page 17. 


